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exploratory  data  analysis  shows  that  low  (high)  dew  point  de- 
pression is  associated  with  the  existence  (non-existence)  of 
low-level  stratus  at  Moffett  Field.   Procedures  for  and  results 
of  various  methods  of  fitting  logistic  models  to  the  data  are 
described.   The  fitted  models  were  used  to  forecast  stratus  op 
reserved  data  sets  (cross-validation).   Results  of  the  cross- 
validation  are  given. 


LOW-LEVEL  STRATUS  PREDICTION  USING  BINARY  STATISTICAL 
REGRESSION:   A  PROGRESS  REPORT,  USING  MOFFETT  FIELD  DATA 

Donald  P.  Gaver      Patricia  A.  Jacobs 

Operations  Research  Department 
Naval  Postgraduate  School 

U .   Executive  Summary 

In  this  paper  various  statistical  models  and  techniques 
are  employed  to  forecast  the  existence  ot  low-level  stratus 
conditions.   They  are  illustrated  for  data  at  an  airport 
(Moffett  Field,  Sunnyvale,  California). 

In  Section  2  the  data  set  is  described  and  the  results 
of  a  preliminary  exploratory  data  analysis  are  given.   These 
suggest  that  dew  point  depression  should  be  predictive  of  the 
existence  of  stratus.   Generally,  low  (high)  dew  point  depres- 
sion is  associated  with  the  existence  (non-existence)  of  stra- 
tus.  This  association  is  also  made  evident  by  a  spectral 
analysis  of  hourly  stratus  levels  and  dewpoint  depression 
described  in  Appendix  F. 

The  remainder  of  this  paper  describes  procedures  for  and 
results  of,  fitting  logistic  models  to  the  data  described. 
Validation  of  the  models  are  addressed  as  well.   The  basic 
logistic  model  is 

expjx  3_} 
P{Y  =  1 (explanatory  variable  x}  =  1  +  exp|x  p J 

where   21   is  a  p-vector  (row)  of  explanatory  variables  and   3_ 
is  a  p-vector  (column)  of  coetticients  to  be  determined. 


Appendix  hi  suggests  several  mathematical  justifications  tor 
use  ot  the  logistic  regression  model. 

We  have  used  various  methods  to  tit  logistic  models  tor 
use  as  predictors  on  reserved  data  sets  (cross-validation). 
Our  cross-validation  experiences  are  reported  in  Appendices  a 
through  D.   Appendix  G  contains  the  asymptotic  distribution  ot 
a  threat  score ,  which  is  one  of  the  statistics  we  use  to  compare 
procedures . 

Appendix  A  reports  on  use  ot  the  stepwise  logistic  regres- 
sion procedure  of  the  BMDP  computer  package.   The  procedure 
chooses  variables  to  be  used  in  the  regression  from  a  menu  of 
variables  given  to  it.   The  BMDP  fits  are  then  used  to  predict 
the  occurrence  ot  stratus  tor  independent  data,  i.e.  from  dif- 
ferent years.   we  find  that  the  stepwise  feature  must  be  used 
with  caution;  it  tends  to  overtit,  inducing  variables  which 
greatly  increase  the  standard  error  ot  the  variables  tirst 
included  in  the  regression.   bucli  overt  ittiny  degrades  the 
predictive  powers  ot  the  model. 

Copas  (iybJ)  points  out  that  a  regression  model,  tit  by 
maximum  likelihood  (or  least  squares)  to  one  set  ot  data,  and 
then  used  tor  prediction  on  another  set  ot  data,  nearly  always 
tits  or  predicts  the  new  set  ot  data  less  well  than  it  does  the 
original  set.   This  phenomenon  ot  shrinkage  can  become  more  pro- 
nounced it  the  original  regression  model  is  tit  using  a  step- 
wise procedure,  which  tends  to  overtit.   Appendix  b  describes 
and  investigates  a  procedure  suggested  by  Copas  to  compensate 
for  shrinkage  on  both  regressions  tit  with,  and  without,  stepwise 


procedures.   In  our  application,  particularly  when  predicting 
changes  from  stratus  to  no  stratus,  the  shrinkage  procedure 
appears  to  help.   However,  it  appears  to  do  less  well  in  predict- 
ing changes  from  no  stratus  to  stratus. 

In  Appendix  C  robust  estimation  procedures  tor  logistic 
regression  are  described  and  carried  out  on  the  Moffett  Field 
data.   These  procedures  are  less  vulnerable  than  maximum  like- 
lihood estimates  to  a  few  outlying  data  points  which  may  not 
agree  with  the  model.   For  the  particular  cross  validations 
performed,  predictions  using  models  fit  with  robust  procedures 
were  no  better  than  predictions  made  with  the  models  tit  with 
maximum  likelihood.   The  models  obtained  are,  however,  system- 
atically different  from  their  classical  counterparts. 

In  Appendix  D  we  investigate  the  predictive  use  of  lo- 
gistic regression  models  that  are  progressively  updated  to  em- 
phasize recent  data.   The  suggestion  is  that  models  fit  with 
data  which  are  closer  in  time  to  the  dates  on  which  forecasts 
are  to  be  made  may  be  more  relevant,  owing  to  changing  condi- 
tions not  represented  in  the  model,  than  a  model  which  is  fit 
with  data  of  several  previous  years.   We  found  that  models 
with  updating  often  did  at  least  as  well  as  models  without  an 
updating  feature. 

In  summary,  we  have  found  that   £n(dew  point  depression  +  1 
appears  to  be  a  consistently  useful  predictor  of  the  occurrence 
of  stratus.   Low  (high)  dew  point  depression  is  associated  with 
stratus  (no  stratus).   There  is  no  one  procedure  or  model,  among 
those  tried  to  date,  that  appears  a  clear  winner.   If,  for 


example,  one  procedure  does  well  in  predicting  changes  from  no 
stratus  to  stratus,  it  will  often  do  less  well  in  predicting 
changes  from  stratus  to  no  stratus.   We  found  that  none  of  the 
procedures  did  as  well  predicting  the  occurrence  of  stratus  in 
196^  as  it  aid  in  1961.   This  suggests  tnat  perhaps  1962  is  not 
described  by  the  present  models  as  well  as  is  1961,  Demg  in- 
trinsically quite  different  from  the  previous  years  195«-61. 
Mooels  and  methods  that  represent  year-to-year  differences  will 
come  under  investigation  in  tuture. 

turtner  work,  with  other  models,  and  with  data  from  other 
locations,  will  oe  undertaken  to  shed  light  on  this  important 
prediction  problem. 


1  .   Introduction  and  Overview 

The  purpose  of  this  paper  is  to  exhibit  the  use  of  sta- 
tistical tools  and  procedures  tor  forecasting  the  existence  of 
low-level  stratus  conditions  at  an  airport.   The  existence  of 
low  stratus  (less  than  or  equal  to  1UUU  ft.)  forces  the  use  of 
different  methods  of  traffic  handling  than  is  the  case  when 
higher  stratus  levels  prevail.   A  low  stratus  condition  tends 
to  inhibit  flight  operations,  so  it  is  desirable  to  torecast 
its  occurrence.   Furthermore,  it  is  of  interest  to  torecast 
such  conditions  on  a  "single-station"  oasis,  making  use  ot  me- 
teorological measurements  available  only  at  the  location— e.y. 
airport--in  question,  in  case  useful  supplementary  information 
is  unavailable. 

The  forecast ing  approaches  described  here  are  statisti- 
cal in  nature,  meaning  that  extensive  data  concerning  the  re- 
ported hourly  stratus  level  at  an  airport  (Moffett  Field, 
bunnyvale,  CA ) ,  together  with  certain  other  meteorological 
measurements  or  parameters  recorded  and  reported  at  that  loca- 
tion, were  used  as  raw  material  tor  the  forecasts.   These  data 
were. used  to  estimate  the  probability  of  low  stratus  during  a 
daily  period;  the  latter  probability  was  estimated  using  a 
logistic  regression  model,  a  tool  that  has  been  found  useful  in 
biological  and  medical  statistics,  and  that  has  been  previously 
applied  in  meteorology;  cf  Brelsford  and  Jones  (1967), 
Gilhausen  (1979),  Gabriel  and  Pun  (1979).   In  a  later  section 
we  present  various  derivations  or  justifications  of  such  a 
model.   Alternative  models  are  also  suggested,  and  the 
usefulness  of  these  will  be  investigated  in  future  work. 


The  usetulness  of  the  logistic  (or  any  other)  model  must 
be  judged  by  its  performance.   We  have  chosen  to  proceed  by  (i) 
fitting  a  model  to  data  for  certain  specific  years  (1958-196U), 
and  then  (ii)  comparing  the  model  predictions  to  actual  occur- 
rences tor  a  completely  different  period  (1961,  1962).   Such  a 
procedure  is  termed  cross  validation;  see  Hosteller  and  Tukey 
(1977)  for  good  general  discussion  and  references.   The  results 
ot  our  cross  validation  are  reported  suDseguent ly .   Another  in- 
teresting and  possibly  useful  approach  is  to  construct  and  test 
an  adaptive,  automatically  up-dating  torecasting  model  witn 
characteristics  similar  to  "exponential  smoothing"  or  "Kalman 
filtering."   Results  of  some  simple  updating  procedures  tor 
torecasting  will  also  be  reported. 

Successful  torecasting  with  the  aid  of  a  model  reguires 
that  the  data  inputs  be  relatively  "clean,"  or  in  basic  con- 
formity with  the  model.   Occasionally  occurring  data  points 
that  are  out  of  line  for  any  reason,  called  out  1 iers ,  or 
influential  values,  can  radically  change  the  values  ot  Lhe 
model  parameters  obtained  from  statistical  fitting  principles 
such  as  least  sguares  (not  used  tor  fitting  our  logistic  model) 
or  maximum  likelihood  (wnich  is  used).   To  check  tor  such 
maverick,  possibly  detrimentally  influential,  values  it  is 
possible  to  proceed  in  several  ways.   One  is  to  successively 
remove  each  data  point  (actually  a  vector  of  response  and 
explanatory  variables)  and  re-tit  the  model,  watching  for 
radical  changes  in  fitted  model  parameters.   This  method  has 
been  programmed  (in  APL,  on  the  NFS  IBM  3033  system)  and 
exercised;  its  detect  is  that  at  present  just  one  data  point 


is  removed  at  a  time,  so  if  several  points  are  mavericks  this 
fact  may  be  overlooked.   Clever  ways  of  automat ica lly  diminish- 
ing the  effects  of  maverick  points  have  been  discussed  by 
Preyibon  (1982);  exploration  of  the  applicability  of  such  ideas 
to  the  present  stratus  prediction  problem  is  currently  underway 
The  methods  and  some  results  are  reported  here. 

Another  approach  to  the  identification  of  maverick  data, 
and  to  the  possible  discovery  of  an  appropriate  model,  is  by 
computer  graphics.   We  have  initiated  the  examination  ot  tne 
low-stratus  data  on  a  pioneering  graphics  facility  at  Stanford 
Linear  Accelerator  (SLaC);  see  an  article  in  science ,  Kolata 
(lytf2),  for  general  description.   Tne  SLAC  system  allows  an 
analyst  views  of  various  three-dimensional  space  projections  ot 
multidimensional  data-clouds.   Such  examination  helps  to  reveal 
the  association  between  certain  explanatory  ("independent") 
variables  and  the  response  ("dependent  variable")  ot  interest. 
For  example,  examination  ot  our  stratus  data  indicated  that 
changes  in  the  explanatory  variable  dewpoint  depression  tended 
to  be  reflected  in  changes  of  response,  i.e.  low  stratus  level 
probability.   This  association  has  physical  basis,  and  dewpoint 
depression  had  actually  been  included  in  earlier  exploratory 
logistic  fits  at  the  suggestion  of  W.  Sweet  ot  NLPKb';  its 
incorporation  into  the  model  considerably  improves  predictive 
performance . 


2 .   The  basic  Data  bet 

The  statistical  methods  used  in  this  study  were  applied 

to  data  furnished  by  to.  Sweet  of  NEPkF,  to  whom  we  are  grateful. 

In  summary,  these  data  consist  ot  reported  hourly  determinations 

of  : 

(i)   stratus  level,  reported  to  be  at  discrete  levels  ot 
100  ft.  separation;  possible  recorded  levels  are 
k  x  1UU  ft.,   k  =  1 ,2, . . . ,9,1U , . . . , "999"  (no  visible 
stratus  )  . 

(li)   east-west  wind  velocity,   V   ,  at  surface,  in  miles 
per  hour; 

(iii)   north-south  wind  velocity,   V   ,  at  surtace,  in  miles 
per  hour; 

(iv)   temperature,  at  surtace,  in  degrees   F; 
(v)   dewpoint,  at  surtace,  in  degrees   F; 
all  at  Moftett  Field,  California,  for  the  months  ot  July, 
August,  and  September  of  the  years  19bb-1962;  later  data  are 
also  available,  and  remain  to  be  analyzed.   Although  other 
measurements,  e.g.  ot  pressure,  are  in  principle  available,  they 
were  not  utilized  in  the  present  analysis.   Nor  were  measure- 
ments from  neighboring  locations  in  the  San  Francisco  bay  area. 
2.1.   The  Forecasting  Exercise  Data  set 

The  raw  aata  described  above  were  adopted  to  the  fore- 
casting exercises  as  follows: 

(a)   Forecasts  are  made  of  the  existence  of  stratus 
level  less  than  1000  ft.  (<_  90U  ft.)  on  any  hour  between 
lUtUU  pm  (221)0)  on  day   t  ,  and  6:00  am  (U6UQ)  on  day 
t  +  1 .   If  hourly-reported  stratus  level  ever  fell  to  a 
level  _<  900  ft.  during  such  a  period  beginning  on  day  t, 
it  is  agreed  to  say  that  stratus  existed  on  day   t;  otherwise 


that  no  stratus  existed  on  day  t.  Denote  by  the  binary 
indicator  variable  y  the  existence  (non-existence)  of 
stratus  on  day   t   according  to  the  above  definition.   Thus 

1   if   stratus  exists  on  day   t, 


0   if   no  stratus  exists  on  day   t. 

Call   y    the  response  (or  dependent  variable)  when  forecasting 
for  day   t.   Note  that  the  observed  values  of  response  on  pre- 
vious days  (y  _,  ,  y   ~ ,  .  .  .  )  are  available  as  assistance  when 
forecasting  for  day   t.   The  above  definition  of  meaningful 
stratus  agrees  with  instrument/no  instrument  landing  rules  at 
airports,  and  is  thus  of  operational  significance. 

Candidate  explana  tory  (independent)  variables  are  these: 

(b)  wind  velocities  at  6:00  pm  (1800)  on  day   t  , 
items  (ii)  and  (iii)  above; 

(c)  temperature   (T  )   and  dewpoint   (D  )   at  surface  at 
6  :  00  pm  on  day   t ; 

(d)  dewpoint  depression,   a\  =  T^  -  D.   at  6:00  pm 

t     t     t 

on  day   t  ; 

(e)  hours  of  stratus   ( H   , )   between  2200  on  the  previous 
day   t-1  and  0600  on  the  current  day   t  ; 

(f)  existence/non-existence  of  stratus   ( y .  , , y .  2  ,  .  .  .  ) 

on  previous  days. 
Let   NS    denote  the  number  of  consecutive  days  of  stratus  in 
a  run  of  stratus  days  that  includes  day  t-1,  the  day  on  which 
the  prediction  is  made.   NNS    is  the  number  of  consecutive  days 
of  no-stratus  in  a  run  of  no-stratus  days  that  includes  day 
t-1. 


No 


te  that  because  of  the  way  in  which  the  response   y 


is  defined,  it  is  legitimate  and  of  interest  to  forecast   y   in 

terms  of   T\  ,  D^ ,  A .  ,  V  (t),  etc.   These  latter  quantities  are 
t    t    t    x 

all  available  at  6:U0  pm  for  forecasts  applying  later,  i.e.  from 
10:00  pm  to  6:00  am  on  the  following  day.   Of  course  many  other 
functions  of  the  hourly  observations  are  candidates  for 
explanatory  variable  status. 


}  0 


3 .   Preliminary  Analysis 

Before  proceeding  to  the  fitting  of  specific  models,  a 
subset  of  the  data  has  been  examined  in  terms  of  simple  summa- 
ries.  Since  the  objective  is  to  forecast,  we  have  divided  (con 
ditioned)  the  data  for  the  years  1958,  1959,  1960  into  four 
groups : 

observations  such  that  Y±.  -\  =  U,  Yt  =  ^  > 
observations  such  that  y  _ .  =  0 ,  y  =  1  , 
observations  such  that   y   .  =  1 ,  y   =  0  , 


Group  00 
Group  01 
Group  10 
Group  11 


observations  such  that   y  _   =  1,  y   =  1  , 


and  have  then  computed  summaries  of  the  observed  distributions 
of  certain  candidate  explanatory  variables.   The  argument  is 
that  a  noticeable  separation  of  such  distributions  when  predict- 
ing  y.   from  the  particular  explanatory  data  suggests  that  the 
variable  in  guestion  may  be  useful  in  forecasting. 

Note  that  we  have  explicitly  used  the  known  stratus  state 
of  the  system  at   t-1   as  one  important  variable,  wishing  to 
make  full  use  of  persistence,  and  to  improve  upon  it.   We  are 
especially  interested  in  the  power  of  explanatory  variables  and 
their  combinations  to  correctly  forecast  changes  in  stratus 
conditions,  e.g.  from   y  _ -.  =  0   (no  stratus  on  day   t-1)  to 
y   =  1   (stratus  on  day   t).   Simple  persistence  forecasting, 
which  predicts   y   =  y   ,   will  never  identify  prospective 
changes  . 

Computer  graphic  analysis  carried  out  at  SLAG,  plus 
physical  insight,  suggest  that  dewpoinf  depression,   A  ,   should 
be  an  effective  explanatory  variable.   Another  useful  variable 


!  1 


seems  to  be  the  hours  of  stratus  observed  on  day   t-1,  denoted 
by   H        There  are  limitless  other  plausible  explanatory 
variables,  as  well  as  combinations  and  re-expressions  (trans- 
formations) of  the  latter,  but  here  we  look  at  only  two.   One 
systematic  way  of  uncovering  predictive  combinations  of  explana- 
tory variables  is  by  use  of  some  form  of  principle  component  or 
factor  analysis;  such  work  is  not  reported  here.   It  seems  pos- 
sible that  a  robust  principle  component  analysis  may  be  informa- 
tive (see  Gnanadesikan  (19"7"7),  or  Campbell  (1982)),  for  the 
existence  of  groups  of  maverick-like  data  have  been  reported  in 
the  overall  data  base.   Clustering  procedures  may  also  be  of 
va lue . 

Tables  1  and  2  give  a  few  useful  summaries  of  the  behavior 
of  the  candidate  explanatory  variables   A    and   H,__-i/*  these 
have  been  developed  for  the  years  1958,  1959,  i960.   The  figures 
in  parentheses  are  natural  logarithms  of  their  counterparts. 
The  log  transformation  is  suggested  to  symmetrize  the  sample 
distribution  (histogram  or  Tukey  stem-leaf  plot),  which  often 
tends  to  appear  positively  skewed  for  the  above  data.   The 
medians  and  quartiles  are  used  instead  of  the  ordinary  means  and 
standard  deviations  because  of  the  possible  non-robust/resistant 
properties  of  the  latter  traditional  mpasures. 

We  can  draw  the  following  conclusions  from  Table  1: 

(a)  corresponding  summary  figures  for  dewpoint  depression 
(0,M,U)  are  rather  stable  from  year  to  year. 

(b)  dewpoinc.  depression  (or  its  log)  should  have  prognostic 
power:   roughly  speaking, 
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TABLE  1 


Observed  Distribution  of  Dew  Point  Depression  (A  ) 


Year 


Lower  Quart i le 
(0) 


Med  ia  n      Upper  Quart i le 
( M  )  ( 0  ) 


19  58; 


19  5  9; 


1960; 


1 

■> 

0 

1 

■* 

1 

1.1 

■> 

1 

0 

->- 

0 

1 

> 

0 

1 

+ 

1 

0 

* 

1 

0 

* 

0 

] 

+ 

0 

1 

> 

] 

I) 

> 

1 

0 

+ 

0 

9(2.2) 
6(1.-79) 
6(1. -79) 
10(2.3) 

8(2.08) 
6(1.-^9) 
"7(1.95) 
10(2.3) 

"7(1.95) 

6(1.-79) 
-7(1.95) 
10(2.3) 


9(2.2) 
■7(1.95) 

8(2.08) 
13( 2.56  ) 

9(2.2) 
-7(1.95) 
9(2.2) 
14( 2.64) 

10(2.3) 
8(2.08) 
8(2.08) 

13( 2.56) 


10(2.3) 
9(2.2) 
8(2.08) 

17(2.83) 

11( 2.4) 
9(2.2) 
10(2.3) 
16(2.^-7) 

13(2.56) 

9(2.2) 
11(2.4) 
16(2.-7-7) 


1  ) 


(b-1)   if  stratus  is  present  at  time  (day)   t-1,  and  it 
A    is  relatively  high  (9  or  above),  a  change  to 
no  stratus  is  indicated,  while  if   A    is  rela- 
tively  low  (below  9)  the  stratus  condition  tends 
to  continue;  on  the  other  hand 
(b-2)   if  no  stratus  is  present  at  time  (day)   t-1,  and 
if   A    is  relatively  high  (10  or  above)  the  no 
stratus  condition  tends  to  continue,  while  if   A 
is  relatively  low  (below  10)  changes  to  a  stratus 
condi  t  ion  become  more  frequent. 
These  results  are  physically  plausible,  and  appear  con- 
sistently, if  not  overwhelmingly  strongly,  in  the  present  data. 

Figures  1  and  2  show  box  plots  of  dew  point  depression 
and  £n(dew  point  depression  +  1)  for  the  years  19b8-60  (cf.  Tukey 
andMosteller  (1977)).   Each  of  the  four  plots  in  the  figures 
contain  only  those  points  for  which   y   ,  =  i  ■+  j  =  y    for 
i,j  =  0,1  .   The  top  (bottom)  edge  of  the  box  is  the  upper 
(lower)  quartile  of  the  data  set;  the  symbol  within  the  box  is 
at  the  median;  the  lines  connect  the  mean;  and  the  circles  out- 
side the  boxes  represent  outlying  data  points. 

It  appears  from  the  top  two  plots  in  each  figure  that  dew 
point  depression,   A   ,  may  have  more  prognostic  value  if  there 
is  no  stratus  the  day  before.   If  there  is  no  stratus  the  day 
before,  then  high   A    appears  to  be  associated  with  persistence 
of  no  stratus.   Since  the  box  plots  do  overlap,  it  is  clear  that 
A    will  not  provide  perfect  prediction. 
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An  exploratory  spectral  analysis  of  hourly  £n(stratus 
height)  and  £n(dew  point  depression  +  1)  for  1958  described  in 
Appendix  F  also  suggests  that  high  (low)  dew  point  depression 
is  associated  with  high  (low)  stratus  height. 

In  Table  2  are  corresponding  figures  for  hours  of  stratus 
on  previous  days. 

TABLE  2 

Observed  Distribution  of  Previous  Days1  Hours  of  Stratus 

Year  Lower  Quart ile      Median      Upper  Quart i le 

(  0  )  (  M  )  ( 0 ) 

19  58;    1+0:  2  (  0  .  "M  4(1.4)  8(2.1) 

1  +  1:  6(1.8)  "Ml. 9)  8(2.1) 

19  59;    1+0:  3(1.1)  4(1.4)  4(1.4) 

1+1:  4(1.4)  6(1.8)  8(2.1) 

1960;    1  +  0:         3.0(1.1)  4(1.4)  4(1.4) 

1+1:  4(1.4)  6(1.8)  8(2.1) 

Again  the  figures  in  parentheses  are  logs. 

Again  some  indications  from  the  table  are  of  interest: 

(a)  corresponding  summary  figures  are  rather  stable,  but 
somewhat  less  so  than  for   A   , 

(b)  relatively  low  values  of  previous  days'  hours  of  stratus 
tend  to  be  associated  with  change  to  no-stratus  condi- 
tion, but  the  tendency  is  rather  weak. 

The  tendency  noticed  above  may  possibly  be  accounted  for  by  the 

fact  that  an  underlying  weather  system  is  passing  over  the 

Moffett  area.   Towards  the  end  of  its  sojourn  there  the  hours 

of  resulting  stratus  tend  to  gradually  decrease  to  zero. 

1"? 


box  plots  for  the  number  of  hours  of  stratus  the  day 
before  when  there  is  stratus,  tor  years  1958-60  appear  in  Figure 
3.  Lach  figure  contains  only  those  points  for  which  the  current 
day  has  no  stratus  or  stratus  respectively.  There  appears  to  be 
an  association  between  a  high  number  of  hours  of  stratus  the  day 
before  and  persistence  of  stratus.  The  association  does  not 
appear  strong,  however. 

Although  the  above  sort  of  analysis  is  interesting,  it 

tails  to  incorporate  the  joint--poss  ibly  interactive—effects 

of  several  variables.   Note  that  no  such  analysis  is  reported 

here  for  the  other  possible  explanatory  variables  related  to 

surface  wind,  namely   V    and   V   .   Somewhat  surprisingly, 

1  x  y  r  ■*    i  > 

these  have  been  tound  to  have  secondary  value  for  the  location 
and  years  under  investigation. 
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APPENDIX  A 
Logistic  Fitting  ana  Cross-Validat ion  usiny  the  bMDP  Package 

In  Appendix  F  we  give  several  mathematical  just  it  icat  ions 
tor  use  ot  the  logistic  regression  model.   In  the  present  Appen- 
dix results  are  given  ot  titting  various  logistic  models  to 
available  Motlett  field  data  tor  years  lyt>b-bl);  they  are  cross- 
validated  tor  years  1961  and  !Vb2. 

Here  the  term  model  refers  to  the  basic  logistic 
representation 

PJ  y  =  l  |  explanatory  variables  xl  =  -= — -^ — 2_1 —         (A-l) 
1 J  d  — '  1  +  expj  xg  \ 

where  x      is  a  p-vector  (row)  ot  explanatory  variables,  and   j3_ 
is  a  p-vector  (column)  of  coefficients  to  be  determined.   The 
bMOP  package  performs  the  fitting,  i.e.  determination  ot   8_ 
from  observations,  by  maximum  likelihood  or  a  closely  related 
method.  It  also  furnishes  Student  t-values  tor  assessiny  the 
statistical  significance  of  the  coefficients  determined,  and 
has  a  step-wise  facility,  which  enters  explanatory  variables  in 
accordance  with  their  judged  explanatory  value.   The  above  pro- 
cedure assumes  that  the  model  is  appropriate  tor  the  data,  a 
practice  that  may  be  dangerous  in  observational  studies,  as  has 
been  pointed  out  by  Pregibon  (ly«j),  who  suggests  some  remedies. 
An  examination  ot  remedies  tor  dealing  with  possibly  "ill- 
fitting"  data  by  the  logistic  is  currently  in  progress,  and  will 
be  applied  to  the  Motfett  Field,  and  other,  data. 

In  the  exercises  reported,  we  have  fitted  195b-iybU  data 
by  logistic  models  using  the  variable  selection  feature.   Two 
types  of  tits  are  considered.   In  one  type  we  condition  on  the 

2  0 


previous  day's  stratus  state;  in  other  words   Pn(x.)   means  the 
probability  of  stratus  on  day   t  ,  given  no  stratus  on   t-1   and 
the  influence  of  explanatory  variables  x    ;  P-i(x.)  means  the 
probability  of  stratus,  given  stratus  on  t-1.   In  the  other  type 
we  have  fit  all  the  data  at  once,  using  an  indicator  variable  to 
identify  stratus  -  no  stratus  days. 

The  predictions  made  are  categorical:   i.e.  if  the  calcu- 
lated p-value  exceeds  0.5,  stratus  is  predicted,  while  if  below, 
no  stratus  is  predicted.   We  have  cross-validated  predictions 
against  the  years  1961  and  1962. 

Model  A-l :   Prediction,  given  no  stratus  the  previous  day 
(y.i  =  0).   The  explanatory  variables  selected  are:   a  constant, 
Jin  (A   +1),   V   .   The  fit  is  as  follows  with  standard  errors  of 
the  fitted  parameters  in  parentheses  below: 


XB  =  6.63  -  3.65  ln(A  )  -  0.08^8  V 
-  t  y 


(1.^0)  (0.^41) 


( 0.0495) 


where   A  =■   A   +  1  . 

The  cross  validations  results  for  1961-19b2  ( F  means  Forecast, 

A  means  Actual)  are  below. 

1961-1962 

Fract  ion 
Correct 

.93 
.2~> 


Fraction  Correct  = 


8  +  6 


88  +  6  +  "?  +  16 


-  =  .80 


.'! 


Note  that  simple  persistence  forecasting  ("tomorrow  is 
the  same  as  today")  for  both  1961  and  19b2  gives  a  fraction  of 
correct  forecasts  equal  to  0.83  (  (  88  +  "? )  /(  88+7+16+4  )  )  ,  which  is 
actually  slightly  better  than  the  logistic  forecast.   However,  the 
present  logistic  model  does  correctly  forecast  about  one-quarter 
to  one-third  of  the  changes  from  no  stratus  to  stratus  correctly; 
persistence  will  never  correctly  forecast  a  change. 


Model  A-2;   Prediction  given  stratus  the  previous  day   (y   ,  =  1) 

The  fit  is 

i 

x  p  =  6.12  -  3.34  £n(A  )  +  0.30  H 


The  variables  selected  were   £n(A   +  1)   and   H   ,  . 


(2.07)  (0.940) 


t 


(0.0893) 


t-1 


where   A   =  A   +  1  . 

The  numbers  in  parentheses  beneath  the  coefficients  are  standard 
errors  based  upon  the  assumption  of  a  correct  model,  maximum- 
likelihood  fitted. 

The  cross  validation  results  are  below. 


1961-1962 


A  \ 

F 

0 

1 

Fract  ion 
Correct 

0 

1 

16 

!  4 

6 
29 

.73 

.6^ 

Fraction  Correct   =  0.69  = 


16    +    29 


14+29+16+6 

In  this  case  the  logistic  model  did  as  well  as  persistence  (0.66) 
in  predicting  stratus  and  no  stratus.   Furthermore,  it  predicted 
73%  of  the  changes  correctly. 
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The  results  of  the  validation  for  1961-1962  of  the  two 
fits  in  Exercises  A-l  and  A-2  are  combined  in  the  following  fable 


U    ■+     0 

88 

0    +    1 

6 

1    +    0 

16 

1  +  1 

29 

Predict  ion 
Success       Failure 


16 
6 

14 


Fraction  Correct 
0.9  3 
0.2^ 
0."?3 
0.6"? 


The  threa  t  score  for  predicting  changes  from   0  +  1  is 


C 


0  +  1 


N0  +  l  +  F0  +  0 


(  6+16  +  "7) 


=  0.21 


(A-2) 


where 


CM  ,  =  the  number  of  correct  predictions  of  change  from 
0  ->•  1  , 


N 


0+1 


the  total  actual  number  of  changes  from   0*1, 


F„  ,\  =  the  number  of  incorrect  predictions  of  no  change 
0>U    0  +    0. 

Similarly  the  threat  score  for  predicting  changes  from 

1    +    0    is 


1  +  0 


16 


Nn     n  +  F,     .  16+6+14 

1  +  0       1  +  1 


0.44     . 


(A-3) 


The  threat  score  for  predicting  all  changes 


TT  = 


C    +  C 
0+1    1*0 


Vl  +  Nl  +  0  +  F0+0+Fl+l 
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The  fraction  of  correct  predictions  using  the  logistic 
models  is 

88+6+16+29 


182 


=  0.-76. 


The  fraction  of  correct  predictions  using  persistence  is 


88+^+29+14 
182 


=  0.^6. 


Of  course  persistence  predictions  will  never  be  correct  when  a 
change  takes  place,  while  the  methods  just  presented,  and  others, 
may  actually  do  quite  well  and  seem  worth  the  extra  trouble. 

Model  A-3:   Prediction  based  on  all  data.   The  variables  selected 
are:   a  constant,   £n(A   +  1),  and  H   , .   The  fit  is 


x  3.  =  6.73  -  3.39  £n(A   +  1)  +   .225  H 


(1.30)   (0.5^0) 


(  0.050"7) 


where   A   =  A   +  1  . 

Again  numbers  in  parenthesis  are  standard  errors. 
The  cross-validation  results  are  below: 


1961-1962 
Pred  ict ion 
Success       Failure 
0  -►  0        91  4 

0+1         5  1*7 

0  ■>  0        16  6 

1+1        30  13 


Fraction  Correct 
0.9  6 
0.2  3 
0.^3 
0.-7  0 
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Fraction    Correct    =    * l    +    5,*    16    +    "°      =    U.7b 

1  o  Z 


traction    Correct    Persistence    =   : — =    u.7b 

lb2 


The  threat  scores  for  predicting  changes  are 


=    U.19 


U  5+17+4 


Tl    =       16    /66    +    13       =    U'46 


TT    =       16    +     5 =    n    34 

11      16+5+17+6+4+13     U,J4 


The  threat  scores  tor  the  fit  using  all  the  data  are  about 
the  same  as  those  tor  the  separate  fits,  i.e.  those  that  condition 
on  whether  or  not  stratus  existed  on  the  day  before.  The  fraction 
of  correct  predictions  of  stratus  and  no  stratus  is  also  about  the 
same  as  that  for  the  two  separate  fits,  and  tor  prediction  by 
persistence.   We  conclude  that  doing  separate  fits  based  on 
whetner  or  not  there  is  stratus  the  day  before  may  not  De 
profitable;  a  single  logistic  model  may  do  as  well  as  two. 
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APPENDIX  b 
bhrinkaye 

The  term  "shrinkage"  is  used  in  connection  with  the 
tollowing  phenomenon:   a  regression  model  tit  by  maximum  likeli- 
hood (or  least  squares)  to  one  set  of  data  which  is  then  used 
for  prediction  on  another  set  of  data  nearly  always  fits  the 
new  set  of  data  less  well  tnan  it  does  the  original  set.   Copas 
(1982)  points  out  that  shrinkage  can  be  more  pronounced  if  the 
original  regression  fit  is  made  with  the  aid  of  a  stepwise  pro- 
cedure;  the  latter  tends  to  overfit.   he  suggests  using  the 
tollowing  logistic  model  tor  binary  prediction: 

P{Y  =  1  j  explanatory  variable  x_] 

n 

exp{6  ,'+K  )       $  !  (x.  -  x  .  )  } 
U     ^  ,    l    l     l 

=   - - (B-l) 

1  +  exp{6 "+K  V   b! (x.-  x. )} 

u   . u ,    1   1    i 
i  =  l 

where   x.   is  the  mean  of  the   i     explanatory  variable  tor  the 
original  data.   {£.}   are  the  MLE  estimators  tor  the  original 
data  and   K   is  a  shrinkage  parameter;   K  =  1  means  that  there 
is  no  shrinkage.   Data-derived  prescriptions  can  be  found  tor 
K,  but  in  the  exploratory  work  reported  here  we  have  found  several 
numerical  trial  values  and  taken  note  of  their  general  effects. 

The  stepwise  regression  procedure  of  BMDP  was  used  to 
fit  a  logistic  model  to  data  t  rom  1958-61).   This  model  with,  and 
without  shrinkage  was  then  used  to  predict  the  occurrence  of 
stratus  in  the  years  1961-62.   Tables  J  and  4  give  the  results 
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of  the  cross  validation.   Note  that  shrinkage  slightly  improves 
the  prediction  of  no  stratus  on  the  following  day. 

Tables  5  and  6  give  the  results  of  fitting  logistic  models 
to  the  data  from  1958-60  and  using  the  models  with  and  without 
shrinkage  to  predict  stratus  in  1961.   Four  different  models  were 
used.   The  parameters  of  the  models  are  as  follows 

Parameters 


Mode  1 


A 

AW 
B 
BW 


constant,  in    A  ,  y  _ -. 

constant,  Jin  A,.,  yt  ,  ,  V  ,  V 

t  J  t  - 1    x    y 

constant,  in    A,  ,  NS.  ,  NNS,  ,  H 


t-1 


constant,  in    A  ,  NS  ,  NNS  ,  H    ,  V  ,  V 


where   A    is  the  dew  point  depression  plus  1. 
The  models  were  fit  using  maximum  likelihood.   Stratus  was 
predicted  on  day   t   if  the  forecast  probability  of  stratus  was 
greater  than  or  egua 1  to   a  .   The  cutoff  point   a   was  taken  to 
be  0.5,  or  alternatively  0.41,  the  fraction  of  days  of  stratus 
during  the  years  1958-1960. 

Tables  7  and  8  give  similar  results  for  the  models  fit  to 
data  in  1958-61  and  validated  on  1962  data.  The  cutoff  point  a 
was  taken  to  be  0.5,  or  alternatively  0.37,  the  fraction  of  days 
of  stratus  during  the  years  1958-61. 

Tables  9  and  10  give  the  threat  scores  for  the  prediction 
of  changes  (equations   A-2,  A-3,  and  A-4 ) . 

The  simplest  model  A  with  a  cutoff  of  0.5  seemed  to  do  as 
well  as  any  of  the  more  complicated  models.   The  rise  of  the 


■I'1 


historical  fraction  of  stratus  days  sometimes  improved  prediction 
of  changes,  but  not  in  all  cases.   The  use  of  shrinkaye  once  again 
often  seemed  to  improve  prediction  ot  changes  from  stratus  to  no 
stratus  but  again  not  uniformly.   Models  A   and  B   with  no 
shrinkage  did  as  well  as  the  stepwise  BMDP  procedure  with  no 
shrinkage . 
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Table  3 
Validation  on  61-62  ot  BhDP  stepwise  Fit 
Using  all  Data  b8-60  with  Shrinkage 


K 

1 

0.6 

!  u. 

,  b 

I  0, 

,4 

""I 

trans  it  ions 

U  -  0 

S 

F 

FC 

S   F 

FC   S 

F 

FC  1  S 

F 

FC 

91 

4 

.96 

93   2 

.9bl 93 

2 

| 
.98 

9b 

0 

1.0 

U  +  1 

5 

17 

.23 

3  19 

.14 

3 

19 

.14 

1 

21 

•ob 

1  ■>  0 

16 

6 

.73 

17   5 

.77 

17 

b 

.77  17 

b 

.77 

1  *  1 

30 

13 

.70 

2b  17 

.60 

26 

17 

1 .60| 22 

21 

.51 

Validation  on  1961  ot  bMDP  stepwise  Fit 
Using  all  Data  58-60  with  Shrinkage 


K 

1 

0 

.6 

0. 

b 

0. 

4 

1 

trans  it  ions 

0  +  0 

s 

F 

FC 

a 

F 

FC 

s 

F 

KC 

S 

F 

FC 

53 

3 

.95 

b4 

2 

.96 

54 

2 

.96 

56 

0 

1   j 

0*1 

3 

8 

.27 

3 

8 

.27 

3 

8 

.27 

1 

10 

.09 

■1  ♦  U 

8 

3 

.73 

9 

2 

.82 

9 

2 

.82 

9 

2 

.  82 

1  +  1 

9 

4 

.69 

8 

5 

.62 

8 

5 

.  62 

7 

6 

.54 

FC  =   traction  correct  predictions 
S  =   number  of  successful  predictions 
F  =   number  of  unsuccessful  predictions 
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Table  4 

Validation  on  iy62  of  bMDP  Stepwise  Fit 

Using  all  Data  38-bl  with  Shrinkaye 


K 

1 

U  .b 

U.5 

0. 

4 

trans  it  ions 
0  +  U 

s 

F 

FC 

S   F 

FC 

b   F 

FC 

s 

F 

FC 

38 

1 

.97 

39    U 

1 

39   u 

1 

8y 

u 

1 

U  +  1 

2 

y 

.18 

U   11 

U 

U   11 

u 

u 

11 

u 

1  -»•  u 

8 

3 

.7j 

8    8 

.78 

a   8 

.73 

8 

3 

.73 

1  ■»  1 

21 

y 

.70 

18   12 

.  bU 

lb   14 

.  88 

11 

ly 

.87 

Model  has  explanatory  variables    constant 


Fst  coetticients 
( Std .  Frrors  ) 


b  .71 
(1.12) 


An(A  ) 

-3.42 
(  U.48y  ) 


H 
t-1 

U  .242 
(  U.U434  ) 


Validation  on  iyb2  ot  Separate  BMDP  stepwise  Fits 
For  Data  Points  with  Stratus  or  No  Stratus  the  Day 
Before  Usiny  data  of  1958-61  with  Shrinkaye 


K 

1 

u 

.b 

U.5 

u. 

4 

trans  it  ions 

s 

F 

FC 

s 

F 

FC 

S   F 

FC 

S 

F 

FC 

U  +  U 

3ti 

1 

.97 

39 

U 

1 

3y   u 

1 

3y 

U 

1 

U  •»  1 

2 

y 

.18 

U 

11 

U 

U    11 

U 

u 

11 

U 

1  ■*  u 

8 

3 

.73 

b 

3 

.  8b 

3     b 

.48 

8 

i 

.73 

1  ->  1 

21 

y 

.711 

22 

8 

.73 

Z4           b 

.  8U 

2b 

•i 

.  8b 

Explanatory  variables  ior  moaei  with  no  stratus  the  day  before 

constant       £n  (  A .  )        V 


Fst  Coetticients 
( S td .  trror  ) 


b  .U3 
( U.42) 


in  (A  ) 

-3.43 
( U.blU ) 


y 

-U.U81 
( U.U4b ) 


Explanatory  variables  tor  model  with  stratus  the  day  betore. 

constant       £n ( A  )        h 

t  t-1 

tst  Coetticients        b . 7 1        -3.3b         U . 2 y 1 

(Std.  Error)         (1.82)        (U.8iy)       (U.U82y) 
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Table  b 
Validation  on  1961  of  Predictions  with  Shrinkage 
of  Models  Fit  with  MLE  Usiny  all  Data  from  1958-6U 


Shrinkage 

Model 

A 

Model 

AW 

Parameter 

Cutoff        0.5 

0.41 

0.5 

0.41 

A\F 

1    0 

FC 

1    0 

FC 

1    0 

FC 

1    0 

FC 

1 

U 

14   1U 
7   60 

.58 
.90 

17    7 
13   54 

.71 
.81 

11   13 
10   57 

.46 
.85 

17    7 
12   55 

.71 
.82 

K  =  1 

trans  it  ions 

S    F 

S    F 

S     F 

S    F 

U  ->•  0 

0  +  i 
i  -»-  0 

1  +  1 

53    3 
3    8 
7    4 

11    2 

.95 
.27 
.64 
.85 

48    8 

5  6 

6  5 
12     1 

.  8b 
.45 
.  54 
.92 

50    6 
3    8 

7  4 

8  5 

1    0 

.89 
.27 
.64 
.  62 

48    8 
5    6 
7    4 

12    1 

.  86 

.45 

.64 

=  s=92 

A\F 

1    U 

1    U 

1 

U 

11   13 
4   63 

.46 
.94 

14   1U 
7   bU 

.58 
.89 

10   14 
4   63 

.42 
.94 

15    9 
11   56 

.  63 
.84 

K  =  U.b 

trans  it  ions 

S    F 

S    F 

S     F 

0  +  0 

0  +  1 

1  +  0 
1  +  1 

54     2 
3    8 
9     2 
8     5 

.96 
.27 
.82 
.  62 

53    3 
3    8 
7    4 

11     2 

.95 
.27 
.64 
.85 

54    2 
3    8 
9    2 
7    6 

.96 
.27 
.  82 
.54 

49    7 
5    6 
7    4 

10    3 

.88 
.45 
.64 
.77 

A\F 

1     0 

1     0 

1     0 

1    0 

1 
U 

8   16 
4   63 

.33 

.94 

14    1U 
7   bU 

.  58 
.89 

8   16 
4   63 

.33 

.94 

12   12 

11   56 

.50 

.84 

K  =  0.5 

trans  it  ions 

S    F 

0  +  U 
0+1 

1  +  U 
.1  +  1 

54    2 
3    8 
y    2 
5    8 

.96 
.27 
.82 
.38 

53    3 
3    8 
7    4 

11    2 

.95 
.27 
.64 
.85 

54    2 
3    8 
9    2 
5    8 

.96 

.27 
.82 
.38 

49    7 
4    7 

7  4 

8  5 

.88 
.31 
.64 
.62 

A\F 

1    U 

1    0 

1    0 

1    0 

1 

U 

5   19 
(J   67 

.79 
1 

14   10 
7   bO 

.58 
.90 

6   18 
1   66 

.25 
.99 

11   13 

8   59 

.85 

k  =  0.4 

trans  it  ions 

S    b' 

S    F 

S    F 

S    F 

0  +  0 
U  +  i 

1  +  0 
1  +  1 

56    0 

1   10 

11    U 

4    9 

1 

.U9 
1 
.31 

53    3 
3    8 

7    4 
11     2 

.93 
.27 
.64 
.  85 

56    0 

Z           9 

1U     1 

4     9 

1 

.18 
.91 
.31 

52    4 
3    8 

7  4 

8  5 

.93 
.27 

.64 
.  62 

FC  =  traction  correct  predictions 

,,  ,,    Number  of  uays  of  stratus  in  195b-19b0 
y  4  j^  =  * 

Total  Number  ot  bays  in  1 


95b-lyb0 


Explanatory  variables  in  nodel   A  =  constant,  y._,/  £n(A  ). 


Explanatory  variables  in  model  AW 


constant,  y,  .,  £n(A  ),  V  ,  V 
J t-1        t     x    y 
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Table  6 
Validation  on  1961  of  Predictions  with  shrinkage 
of  Models  Fit  with  MLE  Using  all  Data  from  19b8-60 


Shrinkage 

Mode 

1  B 

Mode. 

L     BW 

Parameter 

Cutoff        O.b 

0.41 

0.5 

0.41 

A\F 

1    0 

FC 

1    0 

FC 

1    0 

FC 

1    U 

FC 

1 
0 

13   11 
6   61 

.54 
.91 

16    8 
12   55 

.67 
.82 

13   11 
8   59 

.  54 
.88 

16    8 
12   55 

.67 
.82 

K  =  1 

transitions 

S    F 

b    F 

S    F 

b    F 

0  +  0 

0  +    1 

1  +  0 
1  +  1 

53    3 
3    8 
8    3 

10    3 

.95 

.27 
.73 
.77 

48    8 
4    7 
7    4 

12    1 

.86 
.36 
.64 
.92 

52    4 
3    8 
7    4 

10    3 

.93 

.27 
.64 
.77 

48    8 
5    6 
7    4 

11    2 

.86 
.45 
.64 
.85 

A\F 

1    0 

1    0 

1    0 

1 

0 

11   13 
4   63 

.46 
.94 

13   11 
9   58 

.54 
.87 

11   13 
4   63 

.46 
.94 

14   10 
11   56 

.58 
.84 

K  =  0.6 

transitions 

S    F 

b     F 

b     F 

0  +  0 

0  +  1 

1  +  0 
1  +  1 

54    2 
3    8 

y    2 
8    5 

.96 

.27 
.  82 
.  62 

51    5 
3    8 
7    4 

12     1 

.91 

.27 
.64 
.92 

54     2 
3    8 
9    2 
8    b] 

.96 
.27 
.  82 
.62 

49    7 
3    8 
7    4 

11    2 

.  88 
.27 
.64 
.8b 

a\f 

1    0 

1     0 

1    0 

1    0 

1 

0 

11   13 

3   64 

.46 
.96 

13   11 
7   60 

.54 
.90 

10   14 
3   64 

.42 
.96 

13   11 
11   56 

.54 

.84 

K  =  0.5 

trans  it  ions 

S    F 

S    F 

b    F 

b     F 

0  +  0 

0  +  1 

1  +  u 

1  +  1 

55     1 
3    8 
y    2 
8    5 

.98 
.27 
.82 
.  62 

53    3 
3    8 
7    4 

10    3 

.95 
.27 

.64 
.77 

55     1 
3    8 
9    2 
7    6 

.98 
.27 
.82 
.54 

49     7 
3    8 

7    4 
10    3 

.88 
.27 
.64 
.77 

A\F 

1    0 

1    0 

1    0 

1    0 

1 

0 

7   17 
2   65 

.29 
.97 

13   11 
6   61 

.54 
.91 

6   18 
2   65 

.25 
.97 

12   12 
8   59 

.50 

.88 

K  =  0.4 

trans  it  ions 

S    F 

S    F 

S    F 

b    F 

0  +  0 
0+1 

1  +  0 

1  +  1 

56    0 
1   10 
y    2 
b    7 

1 

.09 
.  82 
.46 

53    3 
3    8 
8    3 

10    3 

.95 
.27 
.73 
.77 

56    0 
1   10 
9    2 
5    8 

1 

.09 
.  82 
.38 

52    4 
3    8 

7    4 
9    4 

.93 
.27 
.64 
.69 

FC  =  fraction  correct  predictions 
,       Number  of  Days  of  stratus  in  19b8-l960 
Total  Number  of  Days  in  I9b8-1960 


Explanatory  variables  in  Model   « 
Explanatory  variables  in  Model  bw 


constant ,  ns 


t  ' 
t 


IMNb 


t' 

-1' 


\_i'  £n(V 


constant,  Nb  .  h   .,  £n(A  ),  V  ,  V 
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Table  7 
Validation  Using  1962  of  Predictions  Using  Shrinkage 
and  Models  Fit  with  MLF  Using  all  Data  trora  1958-61 


Shrinkage 
Parameter 

Cutoff 

Model   A 
0.5 

0.37 

Mode] 
U.5 

L   AW 
0.37 

A\F 

1    0 

FC 

1    0 

FC 

1    0 

FC 

1    0 

FC 

1 
0 

26   15 

4   46 

.63 
.93 

29   12 

8   42 

.71 
.84 

25   16 

4   46 

.bl 
.92 

29   12 
7   43 

.71 
.86 

K  =  1 

trans  it  ions 

S    F 

S    F 

S    F 

S    F 

0  +  0 

0  -»■  1 

1  +  0 

1  +  1 

3b    1 
2    9 
8    3 

24    6 

.97 
.18 
.73 
.80 

36    3 
3    8 
6    5 

26    4 

.  92 
.27 
.  55 
.87 

38    1 
2    9 
8    3 

23    7 

.97 
.18 
.73 
.77 

37    2 
J    8 
6    5 

26    4 

.95 
.27 
.55 
.87 

a\f 

1    0 

1    0 

1    0 

1    0 

1 

0 

15   26 

2   48 

.37 
.96 

28   13 

6   44 

.  68 
.88 

15   26 

2   48 

.37 
.96 

26   15 

5   45 

.63 
.90 

K  =  0.6 

trans  it  ions 

S    F 

S    F 

S    F 

S    F 

0  +  0 

0  +  1 

1  +  0 
1  +  1 

39    0 
0   11 
9    2 

15   15 

1 

0 
.82 
.50 

38    1 
2    9 
6    5 

26    4 

.97 
.18 
.  55 
.87 

39    0 
0   11 
9    2 

15   15 

1 
0 

.  82 
.50 

38    1 
2    9 
7    4 

24    6 

1    0 

.97 
.18 
.64 
.87 

a\f 

1    0 

1    0 

1    0 

1 

0 

15   26 

2   48 

.37 
.96 

26   15 

4   46 

.63 

.92 

8   33 
1   49 

.  20 
.98 

26   15 
5   45 

.  b3 
.90 

K  =  0.5 

transit  ions 

S    F 

S    F 

S    F 

S    F 

0  +  0 

0  +  1 

.1  +  0 

1  +  1 

39    0 
0   11 
9    2 

15   15 

1 

0 
.  82 
.50 

38    1 
2    9 
8    3 

24    6 

.97 
.18 
.73 
.80 

39    0 
0   11 

10     1 
8   22 

1 
0 

.91 
.27 

38    1 
2    9 

7    4 
24     6 

.97 
.18 
.b4 
.87 

a\f 

1    0 

1    0 

1    0 

1    0 

1 

0 

7   34 
0   50 

.17 
1  .00 

26   15 
4   46 

.63 
.92 

6   35 
0   50 

.15 
1 

25   lb 
4   46 

.bl 
.92 

K  =  0.4 

trans  it  ions 

S    F 

S    F 

S    F 

S    F 

0  +  0 
0+1 

1  +  0 

1  +  1 

39    0 
0   11 

11    0 
7   23 

1 
0 

1 

.23 

38    1 
2    9 
8    3 

24    6 

.97 
.11 
.73 
.  80 

39    0 
0   11 

11     0 
6   24 

1 
0 
1 
.  20 

38    1 
2    9 
8    3 

23    7 

.97 
.18 
.73 

.77 

0.37 


FC  =  fraction  correct 

Number  of  Days  of  stratus  in  1958-19bl 
Number  of  Days  in  1958-1961 


Model  A  explanatory  variables 


constant,  y 


£n(A  ) . 


t-1'  ~t 

Model  AW  explanatory  variables:   constant,  y.  -.  »  £n(A  ),  V  ,  V 

U    -L  I-  A. 


J  3 


Table  8 
Validation  Using  1962  Data  of  Prediction  Using  Shrinkage 
MLE  Using  all  Data  from  1958-61  and  Models  Fit  with 


Shrinkage 

Mode 

L    B 

Mode. 

L   BW 

Parameter 

Cutoff        U.5 

U.37 

0.5 

U.37 

A\F 

1    0 

FC 

1     0 

FC 

1   0 

FC 

1    0 

FC 

1 
U 

23   18 
4   46 

.56 
.92 

27   14 
8   42 

.66 

.84 

23  18 
3  47 

.56 
.94 

27   14 
7   43 

.66 
.86 

K  =  1 

trans  it  ions 

S    F 

S    F 

S   F 

5    F 

0  +  0 

0  +  1 

1  +  0 

1  +  1 

38    1 
2    9 
8    3 

21    9 

.97 

.18 
.73 
.70 

36    3 
3    8 
6    5 

24    6 

.92 
.27 
.55 
.80 

39   0 
2   9 
8   3 

21   9 

1 

.18 
.73 
.70 

37    2 
3    8 
6    5 

24    6 

1    0 

.95 
.27 
.55 
.80 

A\F 

1    0 

1    0 

1   U 

1 
U 

18   23 
3   47 

.44 
.94 

24   17 
6   44 

.59 
.88 

18  23 

3  47 

.44 
.94 

24   17 
6   44 

.59 
.88 

k  =  U.6 

transitions 

S    F 

S    F 

S   F 

S    F 

0  +  0 
0+1 

1  +  0 
1  +  1 

39    0 
0   11 
8    3 

18   12 

1  .00 

0.00 

.73 

.60 

38    1 
2    9 
6    5 

22    8 

.97 
.18 
.55 
.73 

39   0 
0  11 
8   3 

18  12 

1.00 
0.00 
.73 
.60 

38    1 
2    9 
6    5 

22    8 

.97 
.18 
.55 
.73 

A\F 

1    U 

1    0 

1   U 

1    0 

1 
U 

17   24 
3   47 

.41 
.94 

24   17 
6   44 

.59 
.88 

15  26 
3  47 

.37 
.94 

24   17 
5   45 

.59 
.90 

k  =  U.5 

transitions 

S    F 

S    F 

S   F 

S    F 

0  +  0 
0+1 

1  +  0 
1  +  1 

39    0 
U   11 
8    3 

17   13 

1  .00 

0.00 

.73 

.57 

38    1 
2    9 
6    5 

22    8 

.97 
.18 
.55 
.73 

39   0 
0  11 
8   3 

15  15i 

1 
0 

.73 
_._b0_. 

38    1 

2    9 

7    4 

_22 8 

1    0 

.97 
.18 
.64 
.73 

A\F 

1    0 

1    U 

1   0 

1 
0 

10   31 

3   47 

.24 
.94 

24   17 
5   45 

.59 
.90 

7  34 
3  47 

.17 
.94 

23   18 

4   46 

.  56 
.92 

k  =  U.4 

trans  it  ions 

S    F 

S    F 

S   F 

S    F 

0  +  U 

0  +  1 

1  +  0 
1  -v  1 

39    0 
U   11 
8    3 

1U   2U 

1  .00 

0.00 

.73 

.33 

38    1 
2    9 
7    4 

22    8 

.97 
.18 
.64 
.73 

39   0 
0  11 
8   3 
7  23 

1  .U0 
0.00 
.73 
.23 

38    1 
2    9 
8    3 

21    9 

.97 
.18 
.73 
.70 

FC  =  traction  correct  predictions 
y    -,7  .  Number  of  Days  of  stratus  in  1958-1961 


Number  of  Days  in  1958-1961 
Model  B  explanatory  variables:   constant,  NS 
Model  BW  explanatory  variables:   constant,  NS  ,  NNS  ,  H 


t,  NNS  , 


Ht_1,  £n(At)  . 

1#  *n(A  ),  V  ,  I 
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Table  9 

Threat  Scores  for  1961  Validation 

of  Models  fit  with  MLE  on 

data  from  1958-1960 


Model 

A 

AW 

B 

BW 

Cutpoint 

0.5 

0.41 

0.5 

0.41 

0.5 

0.41 

0.5 

0.41 

T 
0 

0.21 

0.26 

0.18 

0.26 

0.21 

0.21 

0.20 

0.26 

K  =  1    T1 

0.54 

0.50 

0.44 

0.58 

0.57 

0.58 

0.50 

0.54 

TT 

0.37 

0.35 

0.30 

0.39 

0.39 

0.35 

0.34 

0.38 

T 
0 

0.23 

0.21 

0.23 

0.28 

0.23 

0.19 

0.23 

0.17 

K  =  U.6  T, 

0.56 

0.54 

0.53 

0.50 

0.56 

0.58 

0.56 

0.54 

TT 

0.41 

0.37 

0.40 

0.38 

0.41 

0.36 

0.41 

0.32 

Tu 

0.23 

0.21 

0.23 

0.22 

0.25 

0.21 

0.25 

0.17 

K  =  U.5  T1 

0.47 

0.54 

0.47 

0.44 

0.56 

0.50 

0.53 

0.50 

TT 

0.38 

0.37 

0.38 

0.32 

0.43 

0.36 

0.41 

0.31 

T 
0 

0.09 

0.21 

0.18 

0.20 

0.09 

0.21 

0.09 

0.20 

K  =  0.4  T1 

0.55 

0.54 

0.50 

0.44 

0.50 

0.57 

0.47 

0.47 

TT 

0.39 

0.37 

0.39 

0.32 

0.34 

0.39 

0.33 

0.33 

Model 


A 

AW 
B 
BW 

0.41 


Explanatory  Variables 


constant 
constant 
constant 
constant 


*nUt)  yt_1 

*n(At)  yt_1 

£n(A  )  NS 

£n(At)  NS 


V 

x 

NNS, 
NNS, 


V 


t-1 
t-1 


V 


No.  days  of  stratus  during  1958-60 
No.  days  during  1958-60 


V 
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Table  1U 

Threat  Scores  for  1962  Validation 
of  Models  fit  with  MLE  on  data 
from  1958-1961 


A 

AW 

B 

BW 

Cutpoint 

0.50 

0.37 

0.50 

0.37 

0.50 

0.37 

0.50 

0.37 

To 

U.17 

0.21 

0.17 

0.23 

U.17 

0.21 

0.1b 

0.23 

K 

=  1    Ta 

0.47 

0.40 

0.44 

0.40 

0.40 

0.35 

0.40 

0.32 

TT 

0.34 

0.31 

0.33 

0.32 

0.31 

0.29 

0.32 

0.30 

T 
0 

0 

0.17 

0 

0.17 

0 

0.17 

0 

0.17 

K 

=  0.6  T1 

0.35 

0.40 

0.35 

0.41 

0.35 

0.32 

0.35 

0.32 

TT 

0.24 

0.30 

0.24 

0.31 

0.24 

0.26 

0.24 

0.26 

T 
0 

0 

0.17 

0 

0.17 

0 

0.17 

0 

0.17 

K 

=  0.5  T 

0.35 

0.47 

0.30 

0.41 

0.33 

0.32 

0.31 

0.37 

TT 

0.24 

0.34 

0.23 

0.31 

0.23 

0.26 

0  .22 

0.29 

T 
0 

0 

0.17 

0 

0.17 

0 

0.17 

0 

0.17 

K 

=  0.4  T1 

U.  J2 

U.47 

0.31 

U.44 

0.26 

0.37 

0.24 

0.40 

TT 

0.24 

0.34 

0.24 

0.33 

0.19 

0.29 

0.18 

0.31 

Model 


A 

AW 
B 
bW 

0.3"7 


Explanatory  Variables 


constant 
constant 
constant 
constant 


*n(At)  yt_1 
*n(At)  yt_1 

£n(A  )  NS     NNS. 


v 


V 


£n(A  )  NS 


NNS 


t-1 
t-1 


V 


No.  days  of  stratus  during  1958-61 
No.  days  during  1958-61 


V 


y 
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APPENDIX  C 
Robust  Estimation  for  binary  Logistic  Regression. 

Maximum  likelihood  estimates  are  susceptible  to  outlying 
data  points:   they  are  unduly  influenced  by  a  few  (exceptional) 
data  points  which  may  not  agree  with  the  assumed  model.   Pregibon 
(1982)  suggests  robust  procedures  which  yield  estimates  that  are 
resistant  to  a  tew  such  exceptional  data  points.   The  procedure 
that  has  been  used  in  this  report  is  as  follows. 

Let  the  deviance  of  point   i   be 


d.  =  -2  (y  in    p   +  (1-y.)  in(l-p.))f  1  =  1 ,  .  . . , N       (C-l) 


where  in  the  logistic  model 


exp{x  6 } 

1  +  exp{x . 6 } 
—  i 


(C-2) 


and 


x.B  =  BM  +  B,x.,  +  30x    +  ...  +  B  x. 
— l—     0     1  ll     2  i2  p  lp 


(C-3) 


x.,   is  the  value  of  the   k — ■   explanatory  variable  for  the 
lk 

i —   data  point,  and   B    is  the  estimate  of   B.   /  the  regression 

k  K 

coefficient  for  the   k —   explanatory  variable,   x  ;  k  =  l,...,p. 
The  problem  of  finding  the  MLE  estimators  turns  out  to 


be  to  solve  for 


,,,...,B    in  the  non-linear  equations 
u       p 


N 

E    Xik 
i  =  l   1K 


iii 


1  +  e 


Hit 


=  U 


(C-4) 


for   k  =  1 ,  . . . , p . 
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One  possible  robust-resistant  (insensitive  to  outliers) 
procedure  is  to  find  estimators   3  ,...,(3    such  that 


I      w(i)x.k 
1  =  1 


*i 

e 

1 

+   e 

where 


i  f 


w  (  i  )  = 


(H/d.) 


1/2 


=  0 


d  .  <  H 

l  — 


otherwise , 


(C-5) 


(C-6) 


.  th 


d.   is  the  deviance  of  the   i —   data  point  and  the  fitted  model 
l 

at  that  point,  from  (C-l). 

A  value  of   H  =  1.35   was  suggested  by  Pregibon  and  used 
for  the  tuning  constant;  if   H  =  »   the  procedure  carries  out 
the  ordinary  MLR  fitting,  while  as   H   decreases  the  effects  of 
extreme  local  deviance  points  have  progressively  less  effect  on  the 
fitted  model.   Notice  that  the  i —  data-determined  weight,  w(i), 
is  made  relatively  small  if   d(i)   is  large.   Thus  data  points 
which  are  not  well  fit  by  the  assumed  model  will  tend  to  receive 
less  weight  than  others  that  are.   The  resistant  estimates,   3_  , 
are  found  by  iteration.   First  the  MLR  estimate  is  found  and  the 
initial  weights  computed.   Then  (C-5)  is  solved  for   {3  (1), 

K 

k  =  l,...,p}   by  a  Newton-Raphson  procedure.   New  weights   w  (1) 

K 

are  computed  from  (C-6).   Then  these  are  entered  in  (C-5),  and  it 
is  solved  for   (3(2),  k  =  l,...,p};  this  process  repeats  until 
the  iterative  estimates  converge. 
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On  each  day  either  stratus  occurs  or  not.   If  stratus 
occurs  on  consecutive  days  then  a  run  of  stratus  days  is  said 
to  occur.   Let   NS    be  the  length  of  the  run  of  stratus  days 
that  includes  day   t-1.   For  example,  NS   =0   if  the  previous 
day  had  no  stratus,  so   y._-i  =  0?  while 

NSfc  =  2   if   y   1  =  1,  y     =  1,  y     =  0  .   Let   NNS    be  the 
length  of  the  last  run  of  no  stratus  days  that  includes  day   t 


-  1 


Table  (11)  gives  the  estimates  for  five  iterations  of  the 
robust  procedure  applied  to  a  model  using  1958-1960  data.   The 
explanatory  variables  are:   constant,   NS  ,   NNS  ,   H   , ,   £n(A  ) 
here   A    is  the  dew  point  depression  plus  1. 

TABLE  11 
Results  of  Iteration  of  Resistant  Procedure 


w 


Number  of 

I tera t  ion 

Constant 

NSfc 

NNSt 

Vi 

fcn(A  ) 

0  (NILE) 

■ 

6.81 

-0.01 

-0.05 

0.21 

-3.34 

1 

9.30 

-0.04 

-0.05 

0.28 

-4.55 

2 

9.98 

-0.0b 

-0.04 

0.30 

-4.88 

3 

10.16 

-0.0  6 

-0.04 

0.31 

-4.9"7 

4 

10.21 

-0.06 

-0.04 

0.31 

-5.00 

5 

10.22 

-0.06 

-0.0  4 

0.31 

-5.00 

Note  that  except  for  the  estimated  value  of   NNS    the 
resistant  procedure  has  made  the  estimates  greater  in  absolute 
value.   Such  sharpening  of  the  expression  is  a  common  occurrence 
when  robust  logistic  procedures  are  utilized. 
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We  fit  this  model  b  robustly  to  iybb-61)  data  and  then  used 
the  fitted  model  to  predict  the  occurrence  of  stratus  with  a  cutoff 
point  of  U.5.   We  also  robustly  fit  model  B  to  195B-61  data  and 
used  it  to  predict  the  occurrence  of  stratus  in  ly62.   Although  the 
estimated  parameters  using  the  robust  procedure  were  different,  the 
results  of  the  cross-validation  were  almost  the  same  as  with  the 
maximum  likelihood  fit  reported  in  Appendix  b.   Results  of  the 
cross-validations  with  models  fit  robustly  appear  in  Table  lb  at 
the  end  ot  Appendix  D. 
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APPENDIX  U 
Logistic  Models  with  Updating 
Despite  best  attempts  to  develop  a  single  model  with  which 
to  predict  stratus  in  any  given  year,  the  resulting  model  may 
sufter  from  lack  ot  timeliness.   The  basic  reason  is  that  simple 
models  fitted  with  data  t rom  one  period  may  well  not  be  entirely 
relevant  to  another,  owing  to  changing  conditions  not  represented 
in  the  model.   une  attractive  procedure  tor  dealing  with  the  lack 
ot  timeliness  issue  is  to  progressively  update  the  model  fit  so 
as  to  incorporate  recent  data,  i.e.  data  representing  conditions 
near  in  time  to  those  to  be  forecast.   This  is  the  philosophy  ot 
the  well-known  Kalman  filter.   In  the  present  context  the  updating 
procedure  has  been  carried  out  completely  straightforwardly,  i.e. 
by  simply  re-computing  estimates  using  recent  data.   Computation- 
ally economical  and  sophisticated  methods  remain  to  be  developed, 
we  report  the  results  ot  an  investigation  ot  updated  model 
tits  to  predict  the  occurrence  ot  stratus.   Three  updating  schemes 
were  tried. 

1.   a  model  was  initially  tit  using  all  data  from  the  previous 
year.   Then  a  forecast  of  the  occurrence  ot  stratus  was  made  using 
the  model  for  the  first  ten  days  ot  the  current  (forecast)  year. 
These  ten  days  were  then  added  to  the  forecasting  data  set,  ana 
the  eldest,  or  initial,  ten  days  ot  data  were  dropped.   The  model 
was  re-tit  using  the  updated  data.   Using  the  new  model,  the  oc- 
currence ot  stratus  the  next  ten  days  of  the  current  year  was 
forecast.   Then  the  second-eldest  ten-days-worth  of  data  were 
dropped,  and  the  newest  ten  days  were  added,  and  the  model  was 
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re-fit,  forecasts  made,  and  so  the  process  was  continued.  This 
may  be  referred  to  as  a  90-day  rolling  forecast  in  steps  of  ten 
days  . 

II.   A  model  was  initially  fit  using  all  data  from  the  previous 
year.   A  forecast  for  the  occurrence  of  stratus  was  made  for  the 
first  day  of  the  current  year.   This  data  point  was  added  to  the 
forecasting  data  set,  and  the  eldest  point  deleted.   The  model  was 
refit  using  the  altered  modeling  data  set.   A  forecast  of  the  occur- 
rence of  stratus  was  made  for  the  next  day  of  the  current  year. 
This  data  point  was  added  to  the  modeling  data  set  and  the  oldest 
point  was  dropped,  and  so  forth.   This  is  a  rolling  forecast  in 
one-day  steps. 

III.   Same  as  II  but  the  initial  modeling  data  set  includes  only 
the  last  45  points  of  the  previous  year. 

Two  different  sets  of  explanatory  variables  were  tried, 
A   and   B   with  and  without  wind  speeds,  where 
A:   constant,   y   ,,  £n(A  ). 

AW:   constant,   y^  ,,  £n(A),  V  (t),  V  (t) 

t-1       t     x       y 

B:   constant,   NS ,  NNS  ,  H    ,     £n(A  ) 
BW:   constant,   NS  ,  NNS  ,  H   ,,  £n(A  ),  V  (t),  V  (t) 

L  L.        l_    -L  l_  A.  y 

as  before;   A   is  the  dew  point  depression  plus  1. 

A  prediction  of  stratus  was  made  if  the  forecasted  proba- 
bility was  greater  than   a  .   In  most  cases   a  =  0.5.   Additionally, 
a   was  sometimes  taken  to  be  the  fraction  of  the  number  of  days 
of  stratus  over  all  years  previous  to  the  current  year. 

The  results  are  summarized  in  Tables  12-14  of  threat  scores 
(T()fT,iTT)   and  fraction  of  correct  predictions  ( FC  )  .   For  com- 
parison purposes  results  are  also  given  for  prediction  without 
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updatiny.   full  tables  ot  the  numbers  ot  correct  and  incorrect 
predictions  can  be  round  in  Tables  lb-lb. 

As  stated  previously,  the  cutoff  point,   a  ,  for  the  up- 
datiny procedures  was  either  U.b,  or  alternatively,  the  historical 
fraction  ot  days  of  stratus.   for  the  simpler  model  A,  the  use  ot 
the  historical  traction  appeared  to  improve  prediction  ot  stratus, 
but  to  worsen  the  prediction  of  no  stratus.   Usiny  robust  estimates 
in  updatiny  procedure  I  yave  the  same  results  as  usiny  the  simpler 
MLE  estimates.   The  more  complicated  model  b  often  (but  not  always) 
improved  predictions  ot  changes.   Addiny  information  about  winds  to 
either  model  A  or  b  never  improved  prediction  much.   Usiny  shrinkage 
with  the  updatiny  procedure  11  once  again  tended  to  improve  pre- 
diction of  changes  from  stratus  to  no  stratus,  but  tended  to  worsen 
prediction  ot  a  change  from  no  stratus  to  stratus.   Updatiny 
procedure  111  often  seemed  to  do  better  in  predictiny  changes  from 
no  stratus  to  stratus  than  updatiny  procedure  11;  however,  it  did 
worse  when  predictiny  chanyes  t rom  stratus  to  no  stratus.   Upaatiny 
procedure  1  always  did  at  least  as  well  as  in  predictiny  chanyes 
from  stratus  to  no  stratus  but  sometimes  not  as  well  as  III  in 
predictiny  chanyes  from  no  stratus  to  stratus.   Model  b  with  an 
updatiny  procedure  otten  did  better  than  Model  A  with  updatiny 
particularly  in  predictiny  chanyes  from  no  stratus  to  stratus.   In 
summary,  models  with  updatiny  sometimes  did  better  than  models  with 
no  updatiny,  but  the  improvement  was  surprisingly  small. 
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Table  12 

Threat  Scores  for  Changes  and  Fraction  of  Predictions 
Correct  for  iy61  Predictions 

Based  on  Models  With  and  Without  Updating 


Model  A 


Data  Used 
to  Fit  Model 

1958- 

-1960 

1960 

AW 
1958-1960 

Updating 

NO 

NO 

I 

II 

III 

NO 

Method 

MLE 

MLE 

MLE 

MLE 

MLE 

MLE 

a 

0.5 

0.41 

0.5 

0.5 

0.41 

0.5 

0.5 

0.5 

0.41 

To 

TT 

0.21 
0.54 
0.37 

0.26 
0.50 
0.35 

0.26 
0.56 
0.40 

0.25 
0.53 
0.39 

0.39 
0.50 
0.44 

0.25 
0.50 
0.38 

0.29 
0.50 
0.39 

0.18 
0.44 
0.30 

0.26 
0.58 
0.39 

FC 


j  0.81|0.78  ||  0.77  ||0.79  |  0  .  80 |  | 0  .  78 |  | 0 . 78 |  | 0  .  75 | 0 . 79 


Model  B 


Data  Used 
to  Fit  Model 

1958-1960 

1960 

BW 

Updating 

NO 

NO 

I 

11 

III 

1958-1960 

Method 

MLE 

Robust 

MLE 

MLE 

Robust 

MLE 

MLE 

MLE 

a 

0.5 

0.41 

0.5 

0.5 

0.5 

0.5 

0.5 

0.5 

0.5 

0.41 

Tl 

TT 

FC 

0.21 
0.5^ 
0.39 
0.81 

0.21 
0.58 
0.35 

U.^8 

0.21 
0.5^ 
0.39 
0.81 

0.35 
0.56 
0.45 
0.80 

0.43 
0.60 
0.52 
0.84 

0.43 
0.60 
0.52 
0.84 

0.29 
0.56 
0.43 
0.81 

0.29 
0.44 
0.36 
0.77 

0.20 
0.50 
0.34 
0.79 

0.26 
0.54 
0.38 
0.78 

Model  A  has  explanatory  variables:   constant,   y   .,  £n(A  ) 

Model  B  has  explanatory  variables:   constant,   NS  ,  NNS  ,  H   -.  ,  £n(A  ) 

Fraction  correct  using  persistence  is   0.76 
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Table  13 

Threat  Scores  for  changes  and  traction  of  predictions  correct  tor  1962 

Predictions 

Model  A 


Data    Used 
to   Fit   Model 

1958- 

-1961 

1961 

AW 
1958-1961 

Updat  iny 

NU 

NO 

1 

II 

III 

NO 

Method 

MLE 

MLE 

MLE 

MLE 

MLE 

MLE 

a 

O.b 

U  .37 

U.5 

O.b 

U.37 

U.b 

U.b 

U  .b 

U  .37 

Tl 

TT 
VC 

U.17 

U  .47 
0.  J4 

u.7y 

U  .21 
U  .40 
I). 31 

U  .78 

U 

U.47 
0.29 
0.78 

U 
0.31 

U  .lb 

U  .74 

U.21 
U.17 

u.iy 

U.77 

U 

U.31 
U.  lb 
U.7b 

U.14 
U  .27 
U  .  21 
U.7b 

U.17 

U  .44 

U.J3 

|U.78 

U  .23 
U  .40 
0  .32 
0  .79 

Model  b 


Data    Used 
to    Fit    Model 

1958-1961 

1961 

bw 

Updat iny 

NO 

NO 

I 

11 

III 

19b8-19bl 

Method 

MLE 

Kobust 

MLE 

MLL 

Robust 

MLE 

MLE 

MLE 

a 

O.b 

0.37 

O.b 

0.5 

O.b 

O.b 

O.b 

O.b 

O.b 

0.37 

To 

Tl 

TT 
FC 

0.17 
0.40 
0  .31 
0.76 

0.21 
0.40 
0.31 
0.7b 

0.17 
0.3b 
0.28 
0.7b 

0.17 
0.29 
0  .24 
0.7b 

0.1b 
0.37 
0  .28 
0.74 

0.1b 
0.37 
0.28 
0.74 

0.08 
0.26 
0.19 
0.71 

0.20 
0  .22 

0  .21 
0.71 

0.18 
0  .40 
0  .32 
0.77 

0.23 
0  .3b 
0.  30 
0.77 

Model  a  has  explanatory  variables:   constant,   y   ,,  lv\    (a  ) 

Model  b  has  explanatory  variables:   constant,   1Mb  ,  NNb  ,  H   .,  £n  (A  ) 

Fraction  correct  usiny  persistence  is   0.7b 
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Table  lb 

Validations  for  Rolling  tits 
(10  days  at  a  time  initiating  with  only  the 
previous  year ) 


Model 
(Fitting  Procedure) 
Cutoff 

O.b 

i 

(Ml 

a 

B 

(mle; 

O.b 

( Robust ) 

O.b 

ybb-by 

o.b 

0.52 

A\F 
Stratus          1 
no  Stratus       0 

1      0 

FC 

1       0 

FC 

1    u 

FC 

1     0 

FC 

12     12 
b    bb 

.50 

11     13 

b     bb 

.46 
.  bb 

14   10 
b   bb 

.  bb 
.  bb 

14    10 
b   bb 

.bb 
.  btt 

0  +  0 
transitions   U  •*■    1 

1  +  U 
1  +  1 

S      F 

S       F 

s 

S     F 

bl      0 
2     12 
7      8 

10      0 

1 
.17 
.47 

1 

bl       0 
2      12 
7      b 
9       1 

1 
.17 
.47 
.90 

bO    1 

b   y 

b    7 

y   l 

.yb 

.  3b 
.b3 
.90 

bO    1 

b   y 

b    7 

y   1 1 

.  9b 
.  3b 
.53 

.yo 

=)by-b0 

0  .b 

0.267 

O.b 

0  .b 

A\F 

1 

0 

1      0 

1       0 

1    0 

1    0 

22    19 

y   40 

.54 
.  b2 

36       b 

22      27 

.bb 
.bb 

20   21 
b   43 

.49 

20   21 
6   43 

.49 
.  88 

0  +  0 

0  +  1 

1  +  0 

1  +  1 

s    f 

S       F 

S    F 

S    F 

31     3 

b    10 

y    6 
17    y 

.yi 

.67 
.60 
.6b 

23     11 

10       b 

4      11 

26       0 

.  6b 

.67 

.27 

1 

31    3 

4   11 

12    3 

16   10 

.yi 

.27 
.  ttO 
.  b2 

31    3 

4   11 

12    3 

16   10 

.91 

.27 
.80 
.  b2 

460-61 

O.b 

0.41 

O.b 

O.b 

A\F 

1 
0 

1     0 

1       0 

1    0 

1    0 

13    11 
b    bb 

.54 

.  b8 

17       7 
11     bb 

.71 
.b3 

lb    y 
b   61 

.  b3 

.y2 

lb   y 

b   61 

.63 
.92 

0  +  0 

0  +  1 

1  +  0 

1  +  1 

S      F 

S       F 

S    F 

S    F 

bO     b 
4     7 
b     3 
y     4 

.91 
.36 
.73 
.69 

4b       7 
7       4 
7       4 

10       3 

.87 

.64 
.64 
.77 

b2    3 
b    b 
y    2 
y   4 

.yb 

.bb 
.  b2 

.by 

b2    3 
6     b 
y    2 
y    4 

.yb 

.bb 
.  82 
.by 

461-62 

O.b 

0  .37 

O.b 

O.b 

A\F 

1 

0 

1     0 

1       0 

1    0 

1     0 

2b    13 
10    3y 

.  6b 

.bo 

32       9 
12      37 

.7b 
.76 

24   17 
b   43 

.by 

24   17 
b   43 

.by 
.  88 

U  +  0 

0  +  1 

1  +  0 

1  +  1 

S      F 

S       F 

S    F 

S     F 

3b      3 
0     11 
4     7 

Zti              2 

.  y2 

0 
.  36 
.93 

^b                3 
3       b 

2     y 
2y     i4 

.y2 

.27 
.lb 

,y7 

3b    2 

z        y 

7    4 
2  2 b 

.yb 
.ib 

.64 
.73 

3b     2 

2   y 

7     4 
2  2_  _b 

.yb 

.18 
.64 

.73 

number  ot  days  ot  stratus  during  all  previous  years 
number  ot  days  in  ail  previous  years 

3del  A  explanatory  variables:   constant,  y   ,,  £n(A  ) 

Ddel  b  explanatory  variables:   constant,  NS  ,  NINS  ,  H   ,,  £n(A  ) 
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Table  lb 

une  year  Validations  lor  Updatiny  MLE  tits 
of  models  for  one  day  ahead  and  dropping 
the  oldest  day  (cutoff  -  U.b) 


fybb-by 


Model 


Initiating 
Data  bet 


lyby-bU 


A\t' 

1 

U 


A 


transitions 


U 

1 
•  u 
_1_ 

a\f 
1 

II 


1 


ii 
a 


b2 

2 
1 

y 


u  ■>  u 

U  >  1 

1  ->  u 

1  +  1 


24 
1U 


32 

b 

a 
iy 


lybu-bi 


A\F 

i 
u 


12 

a 


b 


u  -*■  u 

U  ■»•  1 

1  +  u 

1  +  i 


1961-62 


A\F 

1 

u 


bl 

4 

a 
a 


u 


lb 
by 


28 

y 


b 


b^ 
U 

4 

2a 


u 

12 

a 
_l 

u 


FC 


4b 

aa 


u 


i 

.14 
.47 

.yu 


17 

4U 

F 


2 

1U 

a 


.by 
.au 


lb 
b2 


s 


F 


bU 
2 

12 
6 

1 


2 

12 

3 

4 

U 


FC 


.33 

.y3 


.yb 

.14 

.au 

.bO 


u 


lb 

a 


y 
by 


.y4 

.33 

.bU 
.73 


U 


12 

by 


24 
12 


17 

ba 


31 
4 
7 

2U 


11 

y 

b 


.by 

.7b 


bl 

a 
lu 


FC 


.63 

.  aa 


u 


a 

1U 


lb 
b7 


F 


22 

a 


iy 

42 


U 


bU 

aa 


13 

41 


.yi 

.36 
.73 
.62 


lb 

y 


n 

bb 


.yi 

.27 
.  44 
.77 


.b4 
.87 


bU 

b 

a 

a 


.  aa 
.  a2 


.ay 

.4b 
.7b 
.  b2 


52 

4 

LU 

la 


2 

11 

b 

a 


.98 

.64 

.bb 

1 


.b4 
.84 


47 
4 

1U 
4 


b 

1U 

b 

0 


FC 


.55 
.8b 


24 
11 


17 

by 


s 


u 


12 

b 


12 
62 


.94 
.27 
.  bb 
.69 


.  bU 

.yb 


31 
4 
8 

21) 


.7b 
.94 
.67 
.40 


.b9 
."78 


b 
11 

a 

b 


U 


1U 


11 

b7 


.91 
.27 
.bU 
.77 


.b4 
.8b 


28 
1U 


lb 
4U 


2 
11 

"7 
2 


.9b 

U 

.5b 

.93 


3b 

2 

4 

2b 


ba 
au 


bb 

4 

y 
a 


23 
8 


18 
42 


.9b 
.bb 
.  82 
.  b2 


bU 
b 
7 
8 


.92 
.18 
.5b 
.87 


37 
1 
b 

22 


2 

1U 

b 

a 


.bb 

.84 


2b 
11 


lb 
b9 


.89 
.4b 
.64 
.  b2 


F 


.9b 

•  U9 
.4b 


3b 
3 

4 


. 73^23 ===7 


.63 

.78 


.90 
.27 
.5b 
.77 


F:   entire  previous  year  used  to  tit  initial  model 
H:   halt  previous  year  used  to  tit  initial  model 
Model  A  explanatory  variables:   constant,  yt_j_»  £n(At) 
Model  b  explanatory  variables:   constant,  N^,  NNbt,  Ht_L'  ^r'(At) 
FC  =  traction  correct  predict;  ions 
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Table  17 

Validation  of  Updating  of  Model  b  with  Shrinkage. 
The  Model  was  initially  fit  with  entire  previous 
year  and  one  point  from  new  year  added  and  oldest 
point  dropped  in  each  iteration. 


K 

1 

0.6 

0. 

5 

0. 

4 

A\F 

1 

0 

FC 

1 

0 

FC 

1 

0 

FC 

1 

0 

FC 

1 

12 

12 

.50 

10 

14 

.42 

8 

16 

.33 

5 

19 

.21 

1960-61 

U 

5 

62 

.93 

4 

63 

.94 

2 

65 

.97 

1 

66 

.99 

transitions 

S 

F 

S 

F 

S 

F 

S 

F 

0  +  0 

53 

3 

.95 

54 

2 

.96 

55 

1 

.98 

56 

0 

1 

0+1 

4 

7 

.36 

3 

8 

.27 

2 

9 

.18 

1 

10 

.09 

1  +  0 

9 

2 

.82 

9 

2 

.82 

10 

1 

.91 

10 

1 

.91 

1  +  1 

8 

5 

.62 

7 

6 

.54 

6 

7 

.46 

4 

9 

.31 

A\  F 

1 

0 

1 

0 

1 

0 

1 

0 

1 

23 

18 

.56 

19 

22 

.46 

19 

22 

.46 

15 

26 

.37 

1961-62 

0 

8 

42 

.84 

5 

45 

.90 

4 

46 

.92 

4 

46 

.92 

trans  it  ions 

S 

F 

S 

F 

S 

F 

S 

F 

0  +  0 

37 

2 

.95 

38 

1 

.97 

38 

1 

.97 

38 

1 

.97 

0  +  1 

1 

10 

.09 

0 

11 

0 

0 

11 

0 

0 

11 

0 

1  +  0 

5 

6 

.45 

7 

4 

.64 

8 

3 

.73 

8 

3 

.73 

1  +  1 

22 

8 

.73 

19 

11 

.63 

19 

11 

.63 

15 

15 

.50 

FC  =  fraction  correct 

Model  B  explanatory  variables:   constant,  NS  ,  NIMS  ,  H   ,  ,  £n(A.  ) 
Cutoff  =  0.5  . 
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Table  18 

Validations  tor  MLE  fits  without  updating  based 
on  different  amounts  of  historical  data 


Validation  Yr. 
Historical  data 

1961 
196U       1958- 

SO 

1962 
1961      1958-1 

51 

A\F 
1 
0 
transitions 
0  +  U 

0  -  1 

1  -►  0 

1  -»■  1 

1    0 

FC 

1    0 

FC 

1    0 

FC 

1    0 

FC 

Model 

A 

13   11 
10   57 

.54 
.85 

14   10 
7   60 

.58 
.90 

24   17 
3   4"? 

.59 
.94 

26   15 
4   46 

.  63 

.92 

S    F 

S    F 

S    F 

S    F 

48    8 
5    6 
9    2 
8    5 

.86 
.45 
.82 

62 

53    3 
3    8 
7    4 

11    2 

.95 
.27 
.64 

85 

39    0 
0   11 
8    3 

24    6 

1 

0 

.73 

.80 

38    1 
2    9 
8    3 

24    6 

,9~i 
.18 
.73 
.80 

A\F 
1 
U 
transitions 
U  ■>  U 

0  -»•  1 

1  +  u 

1  *  1 

1    0 

FC 

1    0 

FC 

1    0 

FC 

1    0 

FC 

Model 

B 

14   1U 
8   59 

.58 
.88 

13   11 
6   61 

.54 
.91 

22   19 
6   44 

.54 
.88 

23   18 
*4   46 

.56 
.92 

5    F 

S    F 

S    F 

fa    F 

50    6 
6    5 
9    2 
8    5 

.89 

.55 

.  82 

62 

53    3 
3    8 
8    3 

10    3 

.95 
.27 
."73 
.77 

38    1 
2    9 
6    5 

20   10 

.97 
.18 
.55 
.67 

38    1 

2    9 

+  8    3 

21    9 

.97 
.18 
.73 
.70 

Cutoff  point  =  0.5 


Model   B 
the  same 


fit  robustly  to 
results  as  mLF. 


data  1958-60  and  cross-validated  on  1961  yives 


Model   B   fit  robustly  to  data  1958-61  and  cross-val laated  on  1962  yives 
the  same  results  as  MLF  except  in  the  cases  *  and  +;  tor  *  tne  corres- 
ponding numbers  are  5  ana  45;  for  +  the  corresponding  numbers  are  7  and  4 


50 


APPENDIX  E 

Survival  Models;   Relation  to  the  Logistic  Representation. 
E . 1   Preliminary  Models 

Suppose  a  system  occupies  one  of  two  states  tor  a  varying 
("random")  time  period,  then  switches  to  the  other,  and  back. 
Such  events  occur  at  times   t  =  0,1,2,3,....   Such  is  the  case 
with  the  stratus-no  stratus  fluctuation  that  has  been  studied, 
but  is  also  true  of  many  other  weather-related  events,  rainfall- 
no  rainfall  being  a  prime  example. 

We  discuss  several  traditional  stochastic  models  as  a 
preliminary. 
Model  1 :   Markov  Chain 

Let   Y    denote  the  state  variable  of  the  system  at  time 
t.   Suppose  (here   i,j  =  0,1) 


P{Y  = 


t=JlVi=i}  ■  Pi: 


>  0  ; 


(E-l) 


in  particular,  no  further  past  history  is  useful: 

P{Yt=j|Yt_1=i,Yt_2=a,Yt_3=b *t_t-k...)  -  pl3     CE-2) 

tor  all   i,j   and  all   t  . 

There  is  then  a  long-run  or  steady-state  distribution 
{•n  ,71,}   that  satisfies  balance  equations: 


"o^i  =  "lPio  =  (1_,ro)pio 


(E-3) 


SO 


10 


'01 


PlO  +  Poi 


PlO  +  P01 
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If  such  a  model  truly  described  nature,  i.e.  stratus  level  at 
an  airport,  then   tt  ,   could  be  referred  to  as  the  climatolog  ical 
probability  of  stratus,   (Y  =1),  on  a  day.   Such  a  model  may  be 
fitted  to  data:   one  simply  estimates   Pin'   for  example,  by 
the  fraction  of  changes  from  1  to  0  (stratus  to  no-stratus) 
observed  in  an  observational  period.   The  model  does  not  have 
the  capacity  to  incorporate  physical  parameters  or  explanatory 
variables,  such  as  dewpoint  depression. 

Model  2:      Two-State  Renewal  Process 

Let   S   represent  the  generic  length  of  a  stratus  period, 
i.e.  or  number  of  days  throughout  which  there  is  uninterrupted 
stratus   (Y  =1).   Just  before   S  ,  and  just  after,  there  will  be 
periods  of  one  or  more  no-stratus  days;  let  such  a  generic  pe- 
riod be   C   (C   denotes  "clear");  (throughout  the  period 
Y   =  U).   If   {S.}   is  a  sequence  of  statistically  indpendent 
stratus  periods  from  the  same  distribution,  and   {C  .  }   is  a 
collection  of  corresponding  clear  periods,  then  the  time  history 
of  system  state  appears  as  below: 


•  ©  •  •  • 


•  o  •  o 


++ Ci .  <--s2— >  + c3 


bZ 


In  the  long  run, 

E[S] 

lim  P{Y  =1}  = 


t       E[S]  +  E[C 


Mean  Length  of  Stratus  Period . 

Mean  Length  of  Strat.  +  Mean  Length  of  No-Strat. 

The  above  can  be  called  the  climatological  probability  of  stratus 

on  a  day.   Strictly,  the  two-state  renewal  process  model  stipulates 

that  the  sequence  of  stratus  day  periods   {S . }   is  one  of  inde- 

J 

pendent,  identically  distributed  random  variables,  as  is  the 
sequence  of  clear  day  periods   {C  } ;  the  two  sequences  are  mutu- 
ally independent.   The  Markov  chain  model  is  a  special  case  of 
the  two-state  renewal  model  in  which  stratus  periods,  generically 
S  ,  have  a  geometric  distribution  with  mean   E[SJ  ,  and  the  clear 
periods,   C  ,  have  their,  generally  different,  geometric  distri- 
bution with  mean   E[C]. 

Once  again,  this  model  contains  no  direct  accounting  for 
the  possible  influence  of  explanatory  variables  upon  the  proba- 
bilities of  stratus  state  changes. 

E . 2   The  General  survival  Model 

Suppose  a  forecaster  is  in  action  at  time   t  .   He  easily 
notes  the  current  system  state;  suppose   Y   =  0  ,  i.e.  no 
stratus.   He  wishes  to  predict  the  system  state  at   t  +  1.   A 
believer  in  Model  I   will  act  in  an  actuarial  fashion,  computing 
the  conditional  probability  that  the  same  state  will  prevail 
("survival"  occurs),  given  that  the  current  clear  state  has 
lasted  tor   d   days: 
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P{C>d+l |C>d}  =  e 


-hQ(d+l) 


or 


p{Yt+l=1'Yt=U'Yt-l=° Yt-d+i=0'Yt-d   =1....}  = 


(E-6) 


=  P{C<d+l|C>d} 


(E--7) 


=  1  -  e 


-h  (d+1) 


Similarly,  if  stratus  is  now  present   (Y  =1), 


P{Yt+1=l|Yt=l,Vt_1=l Vt_d+1=l,Yt_d   =0,...} 


=  e 


■h1(d+l) 


(E-8) 


the  quantities   h  (d),   h,(d)   may  be  referred  to  as  the  hazards 
associated  with  the  states  in  question,  for 


-h1 (d+1) 
1  -  e  -  hx(d+l) 


i  t 


h, ( d+1 )    is  small 


is  the  conditional  probability,  or,  picturesquely,  hazard ,  that 
a  stratus  period  of  duration  (  "  lif elength "  or  "age")   d   actually 
"dies",  or  changes  to  a  non-stratus  period  at  age   d+1. 
Similarly  when  a  non-stratus  period  is  in  progress,  the  change 
occurs  with  hazard   h„(d+l). 

A  promising  enterprise  is  now  to  enhance  the  above  forecast 
of  survival,  or  death,  at  age   d+1   by  further  relevant  informa- 
tion about  the  physical  environment  of  the  process.   Under  present 
circumstances,  i.e.  when  forecasting  stratus,  one  might  well  use 
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dewpoint  depression   A    as  well  as  previous  days  of  stratus  (or 
no-stratus).   Other  explanatory  variables  might  well  be  appropri- 
ate, and  can  perhaps  be  identified  from  physical  arguments  aug- 
mented by  graphical  or  other  exploratory  techniques. 

In  order  to  utilize  the  hazard  notion  in  a  regression 
context  it  is  convenient  to  put 

hQ  =  exp{xtaQ}  (E-10) 

where  for  instance  the  vector  of  explanatory  variables  might  be 


and 


xt  =  (lrNStfAt,Ht„lft) 


-0    (  601'  *  *  *  '30p} 


(E-ll) 


(E-12) 


is  the  required  system  of  constants.   A  form  such  as  (E-1U)  can 
never  be  negative,  a  minimal  requirement.   Precisely  the  expres' 
sion  (E-10)  has  been  used  by  Cox  (1972)  for  describing  hazards. 
Actually  Cox's  hazard  is  written  as 


X ( t  )exp{xg} 


(E-13) 


Suppose  observations  are  available  on   n   days:   these 
are  of  the  form 


(yt'Xtl'Xt2'*** ,Xtp)  ' 


where,  as  was  mentioned  earlier,  possibly 


tl 


1,  x    =  ln(Afc),  xfc3  =  Ht_lf  xt4  =  NSt(i.e.  #  days  of 


continuous  stratus)  . 
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Note  that  interactions  and  transformations  can  directly  be 
included;  e.g.  simply  put   xfc5  =  xt3xt2  =  Ht-1  x  ln<At)   to 
represent  an  interaction  term. 

Now  the  likelihood  for  the   3^   vector  is 


n    -hQ(x  )  y      "h0(x  )  1-y 
L(3_  ;y_,x)  =   n  Le        J   Ll-e        J 
t  =  l 


(fcI-14) 


taking  logs,  we  get 


n  -hyix  ) 

MB.)  =  I     tythu(xt)  +  (l-yt)£n[l-e 


=  -  I     ly   explx^J  +  (  l-yt)£n(l-exp{xt3.u})  ]      (E-lb 


and  this  can  be  maximized  by  choice  of   3_  .  ,  a  non-linear 
optimization  task.   The  usual  approach  would  involve  differen- 
tiation with  respect  to   $..   ,  and  solving  the  resulting  non- 
linear system  by  a  variation  of  the  Newton-Raphson  method. 
Package  programs  are  available  tor  such  a  task. 

E . 3   The  Logistic  Model  from  Cox's  Model 

Suppose  a  Cox  model  is  under  consideration  tor  describing 
the  probability  distribution  ot  "age  to  death"  or,  in  the  present 
context,  the  survival  ot  a  stratus  (or  no-stratus)  episode  tor 
another  day.   in  a  simple  torm,  the  probability  ot  survival 
through  t  +  1   in  state   j   (j  =  1,U)   given  that  for  the  past   m 
time  periods  state   j   is  in  effect  is 
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Vm'V  E  p{Yt+i=jlVJ'Vi=^---    Yt-lB-i)"J' WJ'^t-  st} 

-h  (m+l) 
=  e    J       =  exp  [-A  ( t  )exp{x  3.}  ]  .  (E-16) 

Ordinarily   X  (t)   is  thought  of  as  a  deterministic  but 
J 

unknown  function  of   t  ,  i.e.  time  since  start  of  the  process. 
In  an  application  to  stratus  forecasting,  and  to  other  weather 
phenomena,  it  may  be  desirable  to  allow  a  dependence  of  the  basic 
hazard  rate  upon   m  ,  the  duration  or  "age"  of  the  current  epi- 
sode (stratus,  or  non-stratus  as  the  case  may  be):   X  (m).   This 
necessitates  a  specification,  either  parametric  or  non-parametric; 
the  Cox  procedure  in  Cox  [1972]  was  to  estimate   X  .  (m)   non- 
para  metrically. 

In  order  to  associate  the  Cox  model  explicitly  with  the 
logistic,  adopt  the  attitude  that   X.(m)   is  actually  random, 
and  is  independent ly  distributed  from  period  to  period,  with  a 
distribution  characteristic  of  the  state.   In  such  a  case  we  can 
do  no  better  than  to  attempt  to  estimate  the  model 

P  (m,x  )  =  E(exp[-A  (m)exp{x  £}])  ,  (E-17) 

where  the  expectation  operator  E(»)  is  over  the  distribution  of 
the  now-random  hazard.   To  be  quite  specific,  allow   X  .(m)   to 
have  the  Gamma  distribution  for   a  .  ,y  .    >  0  , 


P{Xn(m)  <  x}  =  /   e   J   — hT~\ h^  '         (E-18) 


where   a    and  v         characterize  the  hazard  variability  when 
state   j   is  in  effect.   Now  for  this  distribution  the  expectation 
is  explicitly  in  terms  of  the  Laplace  transform: 
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a  .y  y  . 

3     (a  .y)     J 


P    (m,x.  )    =    /°°exp[-y    exp{x    R }  ]    ,       -j y    dy 

J         L         Q  L  1    V  '  -i  '  J 


J. 


Y3 


a  .  +  exp{.x_4-.§.} 

Y  A 


1 


1  +  1   e^ti 
a  • 
J 


(E-19) 


This  is  the  probability  of  survival  in  state   j   for  one  more 
period  (no  change). 

Now  the  probability  of  a  change  is,  using  the  above 
randomizing  model, 

P{Yt+1*j|Yt=j,xt}  =  1  " 


J 


1    *t£ 
1  +  -  e 

a  . 
3 


and,  in  case   y •  =  1   (mixing  by  an  exponential)  we  find 

-i  M 

a  .   e 
_J 


(E-20) 


P(Yt  +  1^lYt  =  :^t=xt}  = 


-1  ^ 
1  +  a  .   e 


(E-21) 


J 


which  is  precisely  the  logistic  regression  model.   It  is  thus 
clear  that  the  logistic  regression  model  can  arise  from  a  plau- 
sible stochastic  mechanism.   Note  that  the  derivation  presents 
an  alternative  to  the  simple  logistic  model  that  incorporates  one 
more  parameter,  thus  possibly  allowing  for  the  better  represen- 
tation of  a  wider  range  of  binary  response  data  than  by  the 
classical  logistic. 
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E  .4   The  Cox  Survival  Model  with  Stable-Law  Random  Hazard. 

It  is  of  interest  to  investigate  other  ways  of  introducing 
auxiliary  randomness  into  the  Cox  proportional  hazard  survival 
model.   This  process  considered  here  represents  model  parameter 
fluctuation  from  day  to  day  (in  the  present  application)  that  is 
not  covered  by  the  simple  representation 

H*t  +  1  =  J  |  \  =  JfXt  =  xt)  =  exp[-  A  exp{xt3_}]  ;        (E-22) 
instead  the  form  of  the  randomized  model  is  obtained  by  insertiny  a 
term  in  the  hazard: 

P{yt  +  l  =  J  I  Yt  =  J,-t=-t'et}  =  exP["  Xet  exp{>Lte_}j  .     (E-23) 

Now   e    is  not  directly  observable  or  estimable  if,  as 
is  assumed,  only  one  observation  on  a  probabi li ty  depending  on 
each  e         is  available.   Effectively  one  observes  the  marginal 
probability  of   *..  +  i  =  J  '  given   ^  =J   and  values  of  the  explana- 
tory variables   X  : 

P{Yt+l=J I Yt=J'-t=Xt}  =  Ee  (exP[_Aet  exp{xt£}J)  .     (E-24) 

Suppose  now  that  e         obeys  a  positive  stable  law  distri- 
bution (see  feller  (1966),  p.  17U).   In  this  case  the  Laplace 
transform  of   e    is  always  the  form 

-se  Y 

E[e    tJ  =  e  (aS)   ,     U  <  Y  <  1  .  (E-2b) 

Unfortunately,    explicit    formulas    tor    the    density    of       e         are 
generally    not    available;    that    for      y    =    1/2       is    an    exception: 

f       (x;  01,-7)     =      T7~>    e  •  (E-26) 

t  /2n(x/a) 


by 


It  follows  generally  and  directly  from  (E-24)  and  (E-z5) 
that  if   e    is  positive  stable  the  marginal  probability  of  one- 
day  survival  is 

=  exp[-(Aa)T  exp{x_tYB_}]  ,       (fc-27) 

once  again  exactly  a  Cox  model  (i.e.  of  the  form  (E-22))    but  now 
with  the  parameters 


A'  =  Ua) 


1 


(E-28) 


Thus  the  particular  Cox  model  discussed  is  completely  insensi- 
tive to  the  type  of  hazard  randomization  introduced  here.   Notice 
that  the  effects  of  the  explanatory  variables  or  covariates, 
x   ,  as  measured  by  the  magnitudes  of  their  coefficients 
IB  +  y|  ,  y  <  1)  ,  becomes  progressively  smaller  as   y  +  u  ;  the 
latter  "shrinkage"  tendency  is  associated  with  greater  and 
greater  "spread"  of  the   e.   distribution  (here  "spread"  cannot 
be  measured  by  variance,  tor  the  latter  tails  to  exist).   It 
follows  that  the  predictive  (in  terms  of  explanatory  variables) 
power  of  a  Cox  model  could  improve  by  reducing  any  tendency 
towards  hazard  randomization  of  the  type  exhibited,  if  such  is 
poss  ible . 

further  work  on  randomized  Cox  models  yielding  binary  time 
series  will  be  reported  elsewhere. 
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APPENDIX  F 

Spectral  Analysis  of  Hourly  Stratus  Levels  and  Dew-Point 
Depression  for  July-September  1958. 

The  data  for  the  height  of  the  stratus  level  are  hourly 
records,  in  units  of  hundreds  of  feet,  of  the  height  of  the 
stratus  layer.   There  are  2208  such  observations.   The  data  is 
integer  valued  with  a  minimum  of  three  and  a  maximum  of  999; 
1410  of  the  observations  are  999  which  denotes  the  category  of 
no  stratus  (infinite  height);  the  next  largest  observational 
value  is  888,  of  which  there  are  62;  all  the  rest  of  the  obser- 
vations are  less  than  or  equal  to  180. 

Logarithms  of  the  stratus  heights  were  taken  to  reduce 
the  range  of  the  data.   Figure  (4)  shows  the   Jin   (normalized 
periodogram)  of  the  transformed  data;  (cf.  Cox  and  Lewis  (1966) 
pp.  99).   If  the  data  are  uncorrelated  and  stationary  then 
the  values  of  the  normalized  periodogram  will  appear  independent 
and  have  the  unit  exponential  distribution.  The  line  is  at  the 
95%  quantile  for  the  maximum  of  1104  independent  unit  exponentials 
The  largest  peak  occurs  at  91.   Other  peaks  occur  at  1  and  276. 
The  peak  at  91  suggests  that  a  24  hour  cycle  may  be  present;  the 
peak  at  276  suggests  an  eight  hour  cycle.   The  peaks  around  1  may 
be  attributable  to  the  dependence  of  the  data.   A  least  squares 
cyclic  fit  for  the  24  and  eight  hour  cycles  was  next  carried  out. 
The  residuals  from  the  fit  were  then  whitened,  using  an  AR2  proc- 
ess. Figure  (5)  shows  the  log  (normalized  periodogram)  of  the  resi' 
duals  following  the  cyclic  fit  and  AR2  whitening.   There  are  still 
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two  values  of  the  periodogram  above  the  quantile  line  at  91  and 
160.   Figure  (6)  exhibits  the  cumulative  periodogram  of  the 
residuals.   If  the  residuals  were  uncorrelated  and  stationary, 
then  the  cumulative  periodogram  would  have  the  same  distribu- 
tion as  the  order  statistics  of  an  independent  sample  of  1104 
independent  uniform  random  variables.   The  Kolmogorov-Smirnov 
statistic  of  goodness  of  fit  is  1.12  (Theoretical  99%  quantile  is 
1.628)  and  the  Anderson-Darling  statistic  is  1.39  (theoretical  99% 
quantile  is  3.8S"7). 

As  a  result  of  the  above,  we  model  the  logarithm  of 
hourly  stratus  heights  as 


in    L   =  (-1.55)sin(~|)  -  0.322  cos(^|) 


(-0.202)sin(^^)  "  0.300  cos(^)  +  At  ; 
o  o         t 


A   =  0.^50A   1  +  O.U^BA   2  +  E 


m* 


where   E    are  stationary  and  uncorrelated  random  varicbles. 

I 
Figure  ( 7 )  shows  the  residuals  E   . 

A  similar  analysis  was  carried  out  on  £n [dew  point 

depression  +  11  (LDPD).   The  data  range  from  0  to  9.21;  the  values 

have  a  discrete  nature,  but  not  as  noticeably  as  that  of  the 

stratus  levels.   The  Jln-per  iodogram  of  LDPD  is  given  in 

Figure  8.   There  are  visible  peaks  at  92,  186  and  276,  as 

well  as  near  1.   The  peaks  at  92,  186,  and  276  suggest  24  hr, 

12  hr,  and  8  hr  cycles,  respectively.   A  least-squares  cyclic 

fit  was  made,  and  the  residuals  from  the  fit  were  once  again 

whitened  with  an  AH2  process. 
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Figure  9  gives  the  cumulative  per iodogram  of  the  residuals 
with  the  Kolmogorov-Smirnov  and  Anderson-Uarling  statistics.  Our 
model  lor  LDPD  is 


LDPDt  =  0.015  sin(^j-)  +  0.019  cos(^) 


+    0.060    sin(-^j-j)    -    0.031    cos(^y|) 


+    0.087    Sin(il^)    +    0.061    cos(^) 


+    b 


b      =    .802    b  +    .092    b  +    E       . 


A  graph  of  the  residuals   {Ef}   is  presented  in  Figure  10 
Note  the  two  large  residuals. 

Next  the  residuals,   (k  }   °£  tne   ^n   (stratus  neight) 
level  were  regressed  on  ^t^     >    tne  residuals  of   Jin   (dew  point 
depression)  using  a  least-squares  procedure  and  the  robust  bi- 
weight  procedure. 


E*  =  0.0005  +  0.2712  E^ 
(.022)    (.076) 


(Least  Squares) 
(Standard  Errors 


E*  =  0.0073  +  0.084  E^ 


( b  iweight ) 


If  the  two  points  corresponding  to  large  LDPD  residuals 
are  deleted  than  the  following  values  for  regression  coefficients 
are  obtained 


E£  =  0.009  +  .3421  Ed 
(.022)  (.086) 


( Least  Squares ) 
(Standard  Errors) 


d  d 

E   =  0.0083  +  0.1372  E 


( biweight ) 


c-,  ■) 


The  positive  slope  of  the  regression  of  the  residuals 
suggest  that  the  larger  the  dew  point  depression,  the  higher 
the  stratus  level.   Since  the  regression  was  performed  on  the 
residuals  of  both  series  after  detrending  and  whitening  the 
relationship  should  not  be  strongly  influenced  by  non-stationary 
and  dependence  effects  in  the  marginal  series. 
The  small  values  of  the  fitted  slopes  suggest  that  the 
relationship  is  present,  but  not  very  strong.   This  relationship 
together  with  the  box  plots  of  Figure  1  provide  evidence  that 
aewpoint  depression  and  the  existence  of  stratus  are  indeed 
re  lated . 

The  residuals  were  also  examined  for  lagged  relationships, 

i  d 

eg.  a  relationship  between   E    and   ^*._-j  •   No  relationships 

were  evident. 
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APPENDIX  G 
Threat  Scores 
In  this  section  we  discuss  the  asymptotic  distribution  of 
a  threat  score. 

Consider  an  event  that  either  does  or  does  not  occur  on 
day   n,  n  =  1 , . . . ,N  . 
Let 

1     if  event  occurs  on  day   n  ; 


Let 


Let 


U     if  event  does  not  occur  on  day   n  . 


1     it  the  prediction  is  made  that  the  event 
occurs  on  day   n  ; 


U     if  the  prediction  is  made  that  the  event 
does  not  occur  on  day   n  . 


N 

(1-YH1-X  ) 
L  ,      n      n 
n  =  l 


be  the  number  of  correct  predictions  of  the  event  not  occurring; 

N 

s,  =    y    y  x    , 

1     L  ,   n  n 
n  =  i 

be  the  number  of  correct  predictions  of  the  event  occurring; 

F0  ■  I,    (1"VXn  ' 

n  =  l 

be  the  number  of  incorrect  predictions  when  no  event  occurs; 
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F,  =   y   Y  (1-X  )  , 
n-1 


the  number  of  incorrect  predictions  when  the  event  occurs. 

The  threat  score  for  predicting  that  the  event  occurs  is 


T  = 


Sl  +  Fu  +  Fl 


(li-D 


Equations  (A-2),  (a-3),  and  (a-4)  give  threat  scores  for  predic- 
ting changes  from  no  stratus  to  stratus,  changes  from  stratus 
to  no  stratus,  and  all  changes  respectively. 
Note  that 


T  = 


1-S, 


\   L\l     IM   / 


(G-2) 


where   f  (  x  ,  ,  x  ) 


1-x.. 


If  there  is  perfect  prediction,  then   S,.+  S,=  N   and  T=i.   It  b,=U 
then  T=U.   In  the  case  of  predicting  changes  from  no  stratus  to 
stratus,  the  threat  score  would  be   U   if  prediction  of  stratus 
is  done  using  only  persistence. 

Assume   (b()'  t*(,  #  b,f  t' .  )   has  a  multinomial  distribution  with 
parameters   N,  y   ,  y  UA  #  ^n'    Yiu*   Asymptotically  as 

/SU   b'u   Sl   b'l\ 
N  ■*  o>,  I  — — '  — — '  — — '  — —  J     has  a  normal  distribution  with  mean 

(ruu'  Yul'  Yll'  Ylu*  and  covariance   matrix   —  I      where  I    eguals 
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Yuu ( 1  Yuu ' 


YUUYU1 
YUUY11 


YU0Y1U 


YUUYU1 

Yoi(1-Yoi} 


YU1Y11 


YU1Y1U 


YUUY11 
YU1Y11 

Yll(1-Yll) 

Y11Y1U 


YUUYiU 


YU1Y1U 


Y1IY1U 

Yiu(1-Yiu} 


(ct.  Bishop  et  al.  (1975)). 


A  Taylor  expansion  of   1   in  (G-2)    yields 


■r    -lii-  +   J-  l?!>  \  YH     (tl  \ 


+  o 


ma 


XU     YUU/'U     Yll) 


It  tollows  trom  an  application  of  the  multidimensional  6  - 


method  (ct.  Theorem  14.6-2  of  bishop  et  al.)  that  as   N  -»■  °°  ,   T 

Yll 
has  an  asymptotic  normal  distribution  with  mean  ^z anc' 


UU 


1   ^ 
variance   —  a    where 

N 


2  X     Yll  YUU 

a  =  Yii  — -T- 

(i_Yuu) 


1-T 


It   y     is  fixed,  then   a    has  a  maximum  at   y, , 

1 


UU 


at 


which  it  has  the  value 


2<1-^0>  ' 


At   Yn  =  U   ana   y  ±  ±  =  1  -  yU(j  ,  o   =  u  . 


^4 


Another  application  of  the  6-method  shows  that  the  trans- 


formed threat  score   arcsin/T   has  an  asymptotic  normal  distribu 

Tl 
1-T, 


Tll 
tion  with  mean  arcsin/^ and  variance 


00 


1  1 


N  4  (  1  -  y   )  * 
v     Y  00  ' 


Thus,  if   Ynn   is  fixed,  then  the  transformed  threat  score 
arcsin/T   has  a  variance  which  does  not  depend  on   y    .   How- 
ever, both  the  threat  score,   T  ,  and  arcsin/T   can  have  large 
variance  if   Ynn   i-s  cl°se  to   1  (which  will  often  be  the  case) 
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